Modeling of the water level swell during depressurization of the reactor pressure vessel of the boiling water reactor in accidental conditions

Prediction of the medium behavior during depressurization of the Reactor Pressure Vessel (RPV) of the Boiling Water Reactor (BWR) is an important aspect during the design of the RPV and involved safety systems. An unexpected depressurization may be a result of the Loss of Coolant Accident (LOCA) due to a break in the steam or feedwater pipe. Ballesteros [2] performed an in-depth analysis of events in nuclear power plants between 2002 and 2013. He examined 921 events data obtained from four databases and divided them into categories, then into families and then into sub-families. One of the event classifications was according to the maintenance type (periodic, predictive, planned and corrective). The analysis [2] indicated that 47% of the events were related to periodic maintenance and the main affected Rafał BRyk Holger ScHmidt thomas mull ingo GAnzmAnn Oliver HeRBSt

sciENcE aNd tEchNology components were valves. The analysis shows that even such extremely strict standards like these in nuclear industry cannot eliminate human deficiencies. Therefore, nuclear reactors and their safety systems should be designed a priori self-regulating and independent of failure of any other component.
Since during the regular operation the RPV is filled with medium under the saturation conditions, the sudden decrease of pressure results in bubbles formation in the region occupied by liquid. Due to the significant difference of densities between the liquid and steam, the bubbles occurrence leads to the elevation of the two-phase mixture. Fig. 1 shows the single-phase water level and the level swell due to bubbles formation. The presence of liquid at the high level of the RPV may be unwanted by virtue of the operation of some safety systems of BWRs. For instance, functioning of the Emergency Condenser (EC) of the KERENA reactor may be influenced.
KERENA is a medium-capacity BWR developed by Framatome GmbH [6]. It utilizes several innovative passive safety systems that ensure safe depressurization and heat removal in case of an emergency or accidents, including LOCA. Fig. 2 shows a cross section through the containment of KERENA.   the EC according to the fundamental principle of communicating vessels. According to this principle, when the water level in the RPV decreases because of the evaporation of liquid (due to the pressure drop or decay heat emitted by the core), the cold water column in the EC inlet pipe will go down as well, resulting in injection of cold water into the RPV and saturated steam into the EC. The saturated steam, due to the contact with the cold water at the secondary side, would condense in the pipes of the EC leading to heat removal and depressurization of the vessel. However, due to the level swell in the RPV, elevated liquid water may enter the EC tubes, worsening thereby the heat transfer between the medium inside the tubes and the liquid water in the Flooding Pool Vessel (FPV). In order to prevent this, the EC should be installed at the appropriate level. This level may be concluded from prediction of the level swell. In order to examine the performance of the passive systems of KERENA, a dedicated test facility was built at the Components Testing Department of Framatome in Karlstein, Germany. The Integral Test Facility Karlstein (INKA) represents the KERENA containment with a volume scaling of 1:24. The RPV of KERENA is represented by the steam accumulator of the Karlstein large valve test facility (GAP -Großarmaturen-Prüfstand). The vessel has the storage capacity of 1/6 of the KERENA RPV. In order to simulate the decay heat of the reactor core, the vessel is fed with steam by the Benson boiler with maximum power output of 22 MW.
The models presented in this paper are dedicated for simulation of the GAP vessel. Due to the availability of data, comparison of results was performed with reference to other, similar facilities [7,21]. [1] and [21] show the current state of knowledge by presenting a comparison of the reference data against well-known codes i.e. ATHLET and RELAP. The models presented in this paper were also verified against results obtained with these validated codes.
The models were developed with the object-oriented, equationbased Modelica modelling language [10] and OpenModelica Connection Editor environment [11]. In order to compare various approaches for void fraction modeling, several exchangeable function classes were developed. Each class contain different model of void fraction. The functions may be called from the interface of the vessel object i.e. the dialog box with all parameters can be opened and the specific void fraction model may be chosen from the list of available models.

System description and governing equations
The analyzed thermal-hydraulic system is a vertical vessel divided into N volumes. Nl of the volumes are filled at the beginning with liquid water, and during the pressure drop they may contain either single-phase liquid or two phases. N-Nl volumes at the top of the vessel are always occupied by single-phase steam. During the pressure drop (e.g. due to the steam release during LOCA) and resulting bubbles formation, Nl volumes would expand due to the emergence of void. Correspondingly, the N-Nl volumes at the upper part of the vessel would shrink. As the steam release proceeds, the ratio of release mass flow and steam production becomes larger than 1, and therefore the two-phase mixture level decreases again. As the two-phase mixture level goes down, the opposite situation takes place, i.e. the Nl volumes shrink and N-Nl volumes at the top of the vessel expand. Hence, the vessel was discretized with the dynamic grid. Fig. 3 illustrates the nodalization of the vessel.
In the model it was assumed that the pseudo-steady-state occurs, i.e. that large volumetric vapor production makes a small difference in the liquid volume and thus, liquid velocity is close to zero. This assumption results also in the lack of liquid mass flows between volumes, so the change of liquid mass in a single volume is equal to generation of vapor in this volume due to the pressure drop and heat transfer between the fluid and the wall of the vessel. With this presumption, each of the control volume is described by the conservation of mass for the liquid and vapor phases: Gen dp Gen wall with: In (1) -(4) as well as in (5) -(8) ṁ, m, h, ρ, V and q wall stand for mass flow rate, mass, specific enthalpy, density, volume and heat transferred from the wall to the liquid, respectively. Subscripts l and v represent liquid and vapor phases, while Gen,dp and Gen,wall denote generation of vapor due to the pressure decrease in the vessel and heat exchange with the vessel wall, respectively.
The vapor generation due to the pressure drop was calculated from the energy balance formulated for each volume: Differentiation of the (5) gives: , , Substituting (1) and (2) into (6), vapor generation due to the pressure decrease may be obtained: The steam production due to the heat exchange between the fluid and the vessel wall was calculated according to: Regarding the first volume at the bottom of the vessel, an additional term in the balances was taken into account. Since at INKA the steam from the Benson boiler is delivered at the very bottom of the vessel and then, the appropriate amount of liquid is taken out to keep the mass balance during simulation of the decay heat, an additional heat input, and corresponding vapor production and liquid depletion were taken into account in balances for the first volume:  As a boundary condition, the pressure decrease characteristics was given at the top of the vessel. It was also assumed, that the mass flow at the very bottom of the vessel (i.e. entering the first volume) is equal to zero.
Component named prescribedHeatFlow1 serves for modelling of additional heat delivery to the first volume of the vessel (according to (9)). Since [7,21] do not consider any extra heat input, this term was set to zero in the analyzed cases.

Heat transfer
Apart from the pressure drop, a factor which has an influence on the production of steam in the vessel is the heat capacity of the vessel and the heat transfer between the wall of the vessel and liquid. In general, this heat transfer may be written as: In (10) α stands for the heat transfer coefficient, D for the internal diameter of the vessel and hvol represents the height of a single volume. In order to calculate the heat transfer coefficient, the Chen correlation for saturated flow boiling heat transfer was utilized [3]: Chen correlation considers that two-phase flow boiling is governed by the saturated nucleate boiling mechanism and the convection mechanism, and that the contributions of both are additive. In (11) α nb is the heat transfer coefficient of nucleate boiling determined by Forster and Zuber [8]: In (12) λ, c p , ρ, σ, μ and h lv represent thermal conductivity, heat capacity, density, surface tension, dynamic viscosity and latent heat, respectively. Excess temperature ΔT e and pressure Δp e are calculated as: The convective term in Chen correlation was the one proposed by Dittus and Boelter [5]: where Reynolds number Re l was defined as: The Reynolds number factor F in (11) was defined as: where TP stands for two-phase. However, since the ratio is a flow parameter only, Chen determined its value empirically as a function of reciprocal Martinelli parameter. The shaded area in fig. 5 shows the experimental data, while the blue curve is the best-fit result.
In literature several correlation for F as a function of X tt may be found. The most commonly used [13,15,23], and the one used in this work is: The suppression factor S in (11) was defined as a ratio of mean superheat to the wall superheat: Chen suggested to show S as a function of the local two-phase Reynolds number Re TP . Similarly to the F factor, Chen defined a plot with experimental data for the suppression factor S. Fig. 6. Suppression factor S in a function of Reynolds number [4] The S value used in the model was determined according to [15], [19]:

Void fraction
In order to model the level swell, it is necessary to calculate the volumetric content of steam in each volume i.e. the void fraction. In this work three approaches were applied. Each was developed as a separate function class and each can be called from the dialog box of the GAP model ( fig. 4).
The first approach was based on the churn turbulent drift-flux correlation proposed by Wallis [22] and an explicit analytic equation for the average void fraction as a function of the dimensionless superficial vapor velocity [17], [18]. The second one was based on the dimensionless analysis and correlation proposed by Kurbatov [14]. The last approach utilizes empirical correlation proposed by Labuncov [16].

Drift-flux model
Drift-flux correlations are correlations obtained from experimental data which relate the relative speed of vapor to the void fraction of the two-phase mixture.
Distribution parameter C 0 is a parameter introduced by Zuber and Findley [23] in order to account for radial gradients of void without integrating radially. The parameter is defined as the ratio of the average of the product of flux times void fraction to the product of average of flux and of void fraction: According to Wallis [22], the parameter value usually lies between 1.0 and 1.5. The value of 1.5 indicates a very large radial gradient. In some previous studies, C 0 was treated as a free parameter adjusted for each scenario to gain agreement with measured parameters. In this paper the correlation proposed by Kataoka and Ishii was utilized [12]: The definition of the drift-flux may be written by utilizing the distribution parameter as: with j v and j l standing for superficial velocities of vapor and liquid. Using the pseudo-steady-state assumption, the liquid velocity is zero, so the general drift-flux equation is: On the other hand, Wallis [22] proposed the following correlation for the drift-flux in the churn flow: In (26) the bubble rise velocity U ∞ is calculated according to Hamathy [9]:  [4] sciENcE aNd tEchNology Thus,combining (25) and (26): Sheppard and Morris [17] shown that, from the above, the following analytic equation for the average void fraction may be obtained: with dimensionless superficial vapor velocity ψ defined as:

Kurbatov correlation
The model proposed by Kurbatov [14] is dedicated for water systems under the pressure of between 1 and 190 bars. It takes into account densities and viscosities of the liquid and vapor as well as the surface tension. The Kurbatov correlation is:

Labuncov correlation
Labuncov [16] calculates the bubble rise velocity with an additional correction coefficient (28). where: With the above, the following equation was implemented into a Modelica class: The Lubancov model is dedicated for water systems under 1-196 bars. It delivers good results for the Bond number higher than 500: where D b is the diameter of a single bubble.

Results
In order to validate the models, experimental data from literature were utilized. The results were compared against two independent references [7], [21]. Traichel [21] presents experimental data obtained at a small-diameter vessel (0.088 m of internal diameter). Apart from the measurements, Traichel [21] shows her results in the scope of level swell calculation with ATHLET.
The second reference is one of a series of tests performed circa 1980 by General Electric [7]. In this case, the measurements were performed at a larger vessel (0.3048 m diameter). The same data were used by Aumiller [1] in his assessment of RELAP5-3D. Thus, in both cases the Modelica results may be verified not only against experimental data, but also against other codes.

Small-diameter vessel
The diameter of the small vessel was 0.088 m, its height was 4 m and the wall thickness was 0.01 m. Fig. 7 illustrates the vessel. Fig. 7. Small-diameter test facility [20] In the fig. 7 hcol represents water collapse level while hsw is the level of two-phase mixture. The hcol is measurable, as it may be concluded from measured pressure difference between the top and the bottom of the vessel and from the density of liquid under the saturation conditions. Although the level swell hsw cannot be measured directly, the Modelica model can be verified on the basis of measurable hsw1. The hsw1 level is evaluated by measurement of the pressure difference between the upper part of the vessel occupied by steam and sciENcE aNd tEchNology the pressure of water at the height of 1 m. Since during the transient, void emerges in the liquid region, liquid water from below the level of 1 m is transported over this level. Thus, measured pressure difference is higher than during the steady-state at the beginning of the test. On this basis and since the saturated liquid and steam densities are known, the void (and hence the swell) of the region from the bottom of the vessel up to the height of 1 m may be calculated.
Apart from the measurements, Traichel [21] presented the total two-phase mixture level calculated with ATHLET. Therefore, the total level swell obtained with the Modelica model was also compared with another code.
For the purposes of the Modelica model validation, measurements obtained during the test 47 were utilized [21]. Fig. 8 presents the pressure decrease in the vessel during the test and the corresponding pressure applied into the model at the top of the vessel ( fig. 4).

Fig. 8. Measured and modelled pressures -small-diameter vessel
As the pressure decreases according to characteristics presented in the fig. 8, bubbles start to form first in the upper region of the liquid and then, as the decrease proceeds, in the whole volume of liquid. On the other hand, since the water flows out of the vessel, the total amount of water depletes and the collapse level goes down. Similarly, the two-phase mixture level decreases after reaching a peak, after which the production-outflow ratio becomes lower than 1.
Since three different approaches for modelling of the void fraction were utilized, each of them will be presented on a separate graph. Fig. 9 presents the results obtained with the drift-flux model. Apart from measurements and simulation results, measurements uncertainties of 0.1 m were marked at the points of the highest discrepancies between the curves. Although the assumed uncertainty may be considered as a conservative presupposition, the calculated values are always within the range of +/− 0.1 m. Fig. 9 indicates that calculated hsw1 is higher than the measured one. This means that drift-flux model shows the tendency to overestimate the level which leads to the conclusion that the developed model is rather a conservative one. Modelica-ATHLET comparison of the total two-phase mixture levels shows, that ATHLET overestimate the level even more than Modelica -particularly in the initial phase of pressure decrease. Fig. 10 presents comparison of [21] against the results obtained with the Kurbatov correlation.

Fig. 10. Kurbatov-correlation-model results
Although the results in the scope of hsw1 presented in the fig. 10 fit very well the experimental data, the liquid collapse level is higher than the measured one and it goes beyond the marked discrepancy of 0.1 m. Correspondingly, the two-phase mixture level is much lower than the one provided by ATHLET.
In the third case -in which Labuncov correlation was appliedthe void fraction underestimation is even larger. Fig. 11 presents the results for the third approach.

Medium-diameter vessel
The medium-diameter test facility was constructed from a length of 12 inch, schedule 80 pipe [7]. According to SI standard, the outside diameter of the vessel was 0.3238 m, the wall thickness was 0.01748 m and the height of the vessel was 4.27 m. Fig. 12 presents the scheme of the test rig. Fig. 12 shows that in the case of the medium-diameter vessel, pressure differences were measured at several heights. The methodology of concluding the void fraction and level swell from the pressure difference measurements was the same as in the small-diameter ves-sciENcE aNd tEchNology sel. Implementation of more than one measurement points allowed for the estimation of the total swell level.
General Electric [7] performed several tests at this facility. Since Aumiller [1] referred his results to the test number 1004-3, this paper will focus on the same test. Fig. 13 illustrates the pressure drop measured during the test and the corresponding pressure applied to the Modelica model.

Fig. 13. Measured and modelled pressures -medium-diameter vessel
In the scope of the medium-diameter vessel, only the measurements of the two-phase mixture level were available. These, together with the measurements uncertainties and results obtained with RE-LAP5 and Modelica are presented in the fig. 14. Fig. 14 shows, that all of the models provide results within the measurements uncertainties. However, similarly to the case of the small-diameter vessel, Modelica drift-flux model provides the level swell characteristics which the most accurately fits the exact measurements. Comparison of the Modelica drift-flux model against the RELAP5 models shows slightly different characteristics. The RE-LAP5 models deliver higher values of the peak at the beginning of the transients. Then, the level goes sharply down. According to Modelica models, the peak is not so strong. Lack of the measurement data in the initial phase of the test forecloses the verification of results at the beginning of the transients.

Conclusions
In this paper mathematical models for predicting water level swell were described and their implementation into the object-oriented Modelica language was presented. The simulations results were compared against measurements and the results obtained with other codes available in literature. The model based on the differentiated churn-turbulent drift-flux correlation delivered the best results. The model does not require adjustment of any free parameters performed by a safety engineer. Since the model was equipped with a module for additional heat input the decay heat of the nuclear reactor core may be simulated. Future work on the model assumes validation of the model against a large-diameter vessel.