archives of thermodynamics

Vol. 33(2012), No. 1, 107-122 DOI: 10.2478/v10173-012-0005-0

CFD simulation of DEBORA boiling experiments

ROLAND RZEHAK* ECKHARD KREPPER

Helmholtz-Zentrum Dresden-Rossendorf, POB 510119, D-01314 Dresden, Germany

Abstract In this work we investigate the present capabilities of computational fluid dynamics for wall boiling. The computational model used combines the Euler/Euler two-phase flow description with heat flux partitioning. This kind of modeling was previously applied to boiling water under high pressure conditions relevant to nuclear power systems. Similar conditions in terms of the relevant non-dimensional numbers have been realized in the DEBORA tests using dichlorodifluoromethane (R12) as the working fluid. This facilitated measurements of radial profiles for gas volume fraction, gas velocity, bubble size and liquid temperature as well as axial profiles of wall temperature. After reviewing the theoretical and experimental basis of correlations used in the ANSYS CFX model used for the calculations, we give a careful assessment of the necessary recalibrations to describe the DEBORA tests. The basic CFX model is validated by a detailed comparison to the experimental data for two selected test cases. Simulations with a single set of calibrated parameters are found to give reasonable quantitative agreement with the data for several tests within a certain range of conditions and reproduce the observed tendencies correctly. Several model refinements are then presented each of which is designed to improve one of the remaining deviations between simulation and measurements. Specifically we consider a homogeneous MUSIG model for the bubble size, modified bubble forces, a wall function for turbulent boiling flow and a partial slip boundary condition for the liquid phase. Finally, needs for further model developments are identified and promising directions discussed.

Keywords: Subcooled flow boiling; Computational fluid dynamics simulation; Heat flux partitioning; Two-fluid model

* Corresoponding author. E-mail address: r.rzehak@hzdr.de

Nomenclature

a - bubble influence factor

Aw - wall area fraction influenced by bubbles

Cp - specific heat capacity at constant pressure, JK-1 kg"3

dB - bulk bubble diameter, m

dW - bubble detachment diameter, m

D - pipe diameter, m

f - bubble detachment frequency, Hz

H - specific enthalpy of vaporization, J kg-3

J - mass flow rate, kgm-2s-1

k - thermal conductivity, Wm-1K-1

N - nucleation site density, m-3

P - pressure, Pa

Qtot - wall heat flux, Wm-2

Qc - heat flux due to single phase convection, Wm-2

Qq - heat flux due to quenching, Wm-2

Qe - heat flux due to evaporation, Wm-2

s - hydrodynamic wall roughness, m

twait - waiting time, s

- temperature, K

- saturation temperature, K

- liquid subcooling, K

- wall superheat, K

- wall temperature, K

- friction velocity, m s-1

Greek symbols

a - volume fraction

S - viscous length scale, m

ATrefd - reference temperature for correlation of dW, K

ATrefN - reference temperature for correlation of N, K

H - dynamic viscosity, kgm-1s-1

p - density, kgm-3

a - surface tension, Nm-1

Subscripts

G - gas L - liquid

1 Introduction

Tsub T

For engineering calculations, currently the most widely used computational fluid dynamics (CFD) approach to model two-phase flows with significant volume fractions of both phases is the Eulerian two-fluid framework of in-

terpenetrating continua (see e.g. [1]). In this approach, balance equations for mass, momentum and energy are written for each phase, i.e. gas and liquid, separately and weighted by the so-called volume-fraction which represents the ensemble averaged probability of occurrence for each phase at a certain point in time and space. Exchange terms between the phases appear as source/sink terms in the balance equations. These terms consist of analytical or empirical correlations, expressing the interfacial forces, as well as the heat and mass fluxes, as functions of the average flow parameters. Since most of these correlations are highly problem-specific, their range of validity has to be carefully considered and the entire model has to be validated against experiments.

For the case of boiling flows, where heat is transferred into the fluid from a heated wall at such high rates that vapour is generated, additional source terms describing the physics of these processes at the heated wall have to be included. According to Kurul and Podowski [2,3] the imposed total heat flux is decomposed into a sum of contributions due to different mechanisms of heat transfer. The individual components in this heat flux partitioning are then modeled as functions of the unknown wall temperature and other local flow parameters. From the ensuing relation, the local wall temperature corresponding to the total heat flux can be determined.

A CFD wall boiling model implemented in ANSYS CFX [31] following these lines was calibrated and validated by several authors (e.g. [4]) against experimental results of Bartolomej [5]. In these tests, subcooled boiling of water at high pressure flowing upwards in a vertical pipe heated from the outside was investigated and measurements of the axial development of void-fraction, wall temperature and cross-sectionally averaged liquid temperature were provided.

The aim of this work is to investigate the applicability of the CFX models to a different set of experiments, namely the DEBORA tests, in order to demonstrate their general validity and identify specific weak points as well as promising directions for further development. To this end, the previous investigation [4] is taken as a starting point. Changes are introduced rather sparingly and only where this is found indispensable or at least results in a significantly improved model accuracy. The results complement our previous investigation [6].

2 The DEBORA boiling tests

A detailed description of the DEBORA test facility and measurements can be found in [7]. In a vertical heated pipe having an inner diameter of 19.2 mm R12 is heated over a pipe length of 3.5 m. The radial profiles for gas volume fraction and gas velocities at the end of the heated length are measured by means of an optical probe. Furthermore, profiles of bubble size at this position are available. In addition, radial liquid temperature profiles as well as axial profiles of the wall temperature are measured by thermocouples.

Two tests at different conditions (listed in Tab. 1) are selected for the present investigation. These and other test cases have been considered previously [8-11] to assess boiling models implemented in the NEPTUNE code. Saturation properties of liquid and vapour for both pressure levels are given in Tab. 2.

Table 1. System parameters for the selected test cases.

Test DEBORA 1 DEBORA 4

Pressure [MPa] 2.62 1.46

Mass flow rate [kgm~2s~-1] 1996 2030

Wall heat flux [kWm-2] 73.89 76.24

Inlet temperature [°C] 68.52 31.16

Table 2. Material properties of R12 at different pressures.

Pressure condition [MPa] 2.62 1.46

Saturation temperature [K] 360 331

Surface tension [Nm_1] 0.00180 0.00465

Enthalpy of vaporization [Jkg 1 ] 293 xlO3 258 xlO3

Vapor density [kgm~3] 172 85.0

Liquid density [kgm~3] 1.02X103 1.18X103

Liquid Cp [Jkg^K-1] 1.42X103 1.11x10s

Liquid viscosity [kgm~-'-s--1] 89.5x10-® 131xl0~6

Liquid thermal conductivity [Wm_1K_1] 0.0457 0.0558

3 CFD models for subcooled flow boiling

The simulations were made using the Euler/Euler two-fluid and heat flux partitioning models of ANSYS CFX 12.1. Some specific model adjustments made will be detailed later on after summarizing the general settings which were taken similar to other works (e.g. [4]).

Different mechanisms of heat transfer from the hot wall to liquid and gas phases are described by splitting the total heat flux into three contributions due to single phase convection, quenching and evaporation. The heat flux due to single phase convection is described by a heat transfer coefficient calculated from the temperature wall function of Kader [12] as described in [13]. The quenching heat flux is calculated from a one-dimensional transient purely conductive solution following [14]. The evaporation heat flux follows from the amount of vapour generated and is expressed in terms of bubble parameters as

where pG is the vapour density and Hlg is the specific enthalpy of vaporization.

For the bubble detachment size dw, detachment frequency f and nucleation site density N, standard correlations provided by ANSYS CFX are used but some recalibration is found necessary for the present application as will be discussed below.

A change of phase may also take place in the bulk of the fluid. The heat transferred during condensation or evaporation processes is calculated according to the correlation of Ranz and Marshall [15]. For the accompanying changes of bubble size a simple parametrization in terms of the local liquid temperature due to Anglart [16] is used in the reference simulations with suitably recalibrated coefficients. Subsequently, the prediction of bubble size profiles is shown to be greatly improved by applying a homogeneous MUSIG model which in addition to condensation/evaporation also accounts for coalescence/fragmentation.

For momentum exchange between the phases, the Ishii and Zuber drag law [17] was used. Furthermore, a lift force according to Tomiyama [18] and a Favre averaged turbulent dispersion force due to Burns [19] were taken into account. Coefficients determining the strength of these interactions are taken at their standard values for the reference simulations, but will be reconsidered in order to show their influence on the model predictions.

Qe = ^PcHLGd^fN ,

The turbulence of the liquid phase was modelled by a shear stress transport model [20] with the effects induced by the bubbles included as a separate contribution to the effective viscosity [21]. Single phase turbulent wall functions are used for the reference simulations but a modified wall function for boiling flow will be shown to improve the agreement with the velocity profiles.

Material properties were taken as constant at the values given in Tab. 2. Calculations were made on a quasi-2D cylindrical geometry, i.e. a narrow cylindrical sector with symmetry boundary conditions imposed on the side faces. The validity of this simplification has been verified by grid resolution studies and by comparison to a 3D simulation representing a 60o sector of the pipe.

An inlet condition was set at the bottom specifying a typical turbulent pipe flow velocity profile. At the outlet at the top a pressure boundary condition was imposed. On the heated walls, no-slip and free-slip conditions were used for liquid and gas velocities, respectively. Fluxes of mass and energy to each phase derive from the heat flux partitioning discussed above.

4 Results

4.1 Model calibration

Bubble detachment size The bubble size at detachment, dw, depends on the liquid subcooling Tsub = Tsat — Tl as observed by Tolubinsky and Kostanchuk [22] for water at atmospheric pressure. The CFX model is based on a fit of their data to a correlation of the form

dw = drefeATrefd , (2)

where the fit parameters dref and ATrefd depend on other variables having an effect on the detachment size such as flow rate, heat flux and material properties at the system pressure. Adjustment to the Bartholomej tests [5] on water at high pressure by results in different parameter values than needed to match the original data (e.g. [4]).

For the present application to R12 at different pressure levels, ATrefd = 45 K has been left unchanged, but dref is set to 0.24 mm for P = 2.62 MPa and 0.35 mm for P = 1.46 MPa in order to obtain the measured bubble sizes close to the wall.

Bulk bubble size To model condensation of bubbles in the bulk of the subcooled liquid Anglart [16] proposed for the buble diameter size, the simple parametrization:

and suggested values of the coefficients for conditions typical to nuclear energy systems. Again these coefficients have to be adjusted for the present application. Keeping the reference subcoolings Tsub,i = —13.5 K and Tsub,2 = 5 K, the corresponding bubbles sizes are changed to ds,i = 0.035 mm and ds,2 = 0.7 mm for P = 2.62 MPa and ds,i = 0.066 mm and ds,2 = 1.2 mm for P = 1.46 MPa in order to match the measured bubble sizes in the pipe centre.

Nucleation site density Correlations for the nucleation site density are most of the time expressed as power laws depending on the wall superheat Tsup = TW — Tsat, i.e.

but a recent compilation [23] shows that vastly different values are required for the parameters Nref and ATrefN to match different data sets. A likely reason for this fact is that nucleation site density is highly dependent on the microscale topography of the boiling surface, which in turn depends strongly on the processes that were used to finish the surface. These processes are very diverse and in most boiling experiments not specifically controlled. Therefore, these parameters have to be adjusted differently for each boiling surface. In addition they are expected to depend on material properties at the system pressure while dependences on other system variables are likely to be weak owing to the small length scales of the heterogeneous nucleation phenomena.

For the present application direct data on the nucleation site densities are not available, but the measured wall temperatures provide sufficient information to determine the necessary parameter values. We keep p = 1.85 and ATrefN = 10 K from previous work [4] but adjust Nref = 3.0 x 107m"2 for P = 2.62 MPa and Nref = 5.0 x 106 m"2 for P = 2.62 MPa in order to match the measured wall temperatures as shown in Fig. 1 for the DEBORA 1 test case. Clearly the value of Nref = 8 x 105 m"2 used in [4] for the Bartholomej tests [5] here gives unrealistically high wall temperatures.

dB,i{Ts-ub — Tsubt 2) + dB,2{Ts-ub — Tsubti)

sub,2 T sub,l

| 2-0 T3

ra x ra

0.0 I .... I .... I

0 5 10 15

superheat [K]

Figure 1. Wall superheat for different values of the parameter Nref in Eq. (4) for the DEBORA 1 test.

4.2 Reference calculations

Results of calculations using the model with the values of parameters described so far are shown in Figs. 2 and 3 for the test cases DEBORA 1 and 4, respectively. As the comparison of measured (symbols) and calculated (lines) profiles at the end of the heated length shows, overall an agreement between experiment and simulation is found which is comparable to other state of the art investigations in the field. Deviations of up to a factor of two occur for the gas fraction in the core region of the pipe and for the bubble size near the wall. The relative deviations in the gas velocities and the liquid temperature are in the 20 to 30% range. This applies equally to both test cases DEBORA 1 and 4. Several other tests have also been simulated (not shown here due to space limitations) with the same system pressures and mass flow rates. In these calculations similar agreement with the data is found using the same calibration of model parameters as presented here.

4.3 Bulk bubble size

The most obvious deviation between simulation and experiment is the mismatch in the bulk bubble size in the near wall region. To improve the

*xx Eáp

.....//REF=8.0e5

— /NREF=1.5e7 A/REF=3.0e7 /-- NREF=5.0e7 •''- ' " WREF=1.0e8 ----NR£F=3.0e8

superheat [K]

Figure 2. Comparison of measured (symbols) and calculated (lines) profiles for the DEBORA 1 test: reference model.

reference model in this respect, a homogeneous MUSIG approach accounting for both coalescence/fragmentation [24] and evaporation/condensation mechanisms [25] has been implemented. For the present application, the bubble detachment diameter dw from the wall boiling model is imposed as a boundary condition at the wall. Only one gas velocity field but 10 size fractions in the range 0-1 mm (DEBORA 1) and 0-1.5 mm (DEBORA 4) are used. Between the size fractions bubble coalescence is modelled according to Prince and Blanch [26] and bubble fragmentation according to Luo and Svendsen [27]. For both processes rate multipliers were applied and adjusted to match the data. Further work is necessary to avoid this calibration procedure.

Increase of bubble size in the direction away from the wall must at least partly be due to coalescence of bubbles since the liquid is subcooled in a large part of the region over which this increase occurs. As shown in Fig. 4 this is possible with reasonable agreement for both test cases DEBORA 1 and 4 with a single set of values for the adjustable multipliers of coalescence and fragmentation rates. A similar comparison of the other measured quanti-

1 330-

0.002 0.004 0.006 0.008 0.01 radius [m]

« x xExp :-CFX-12

- | Tsat

I»*" —'— —" * 1 X —'—

- 1 » * * * * ,,.>.,.1.?.?."..".;-------- —'—1—'—1—'—

■S c (D C E .33

ffl 0.0004

^ 0.0002

0.0000

i-'—I-1—I

0.002 0.004 0.006 0.008 0.01 radius [m]

t?5 x*tr X *

X M «Exp *

-'- — CFX-12 -'-i-'- —T—

0 0.002 0.004 0.006 0.008 0.01 radius [m]

0.002 0.004 0.006 0.008 0.01 radius [m]

Figure 3. Comparison of measured (symbols) and calculated (lines) profiles for the DEBORA 4 test: reference model.

ties (not shown), however, reveals hardly any differences to the reference calculation at all. Therefore other model variations will be explored in the following.

Figure 4. Comparison of measured and calculated values for DEBORA 1 and 4: parametrized bubble size (dashed lines) and homogeneous MUSIG model (solid lines).

4.4 Bubble forces

The gas fraction profile is determined by the non-drag forces which govern the radial distribution of the bubbles. Thus, the deviations observed in the reference calculations may be related to an improper balance between lift and turbulent dispersion forces. Examining the experimental support for the values of the coefficients that adjust the strength of these forces, it is found that especially the basis for the latter is rather weak. Therefore we have performed additional calculations where the turbulent dispersion force has been changed to the half and double of the reference value.

As shown in Fig. 5 by changing the strength of this interaction, the agreement with the measured gas fraction profile may be improved either in the center of the pipe or near the walls but not at both locations at once. Concerning the other measured quantities, there is a corresponding effect on the velocity profiles while the liquid temperature and consequently also the bubble size are unaffected. These findings suggest that further effort should be directed towards a more precise modelling of the turbulent dispersion force.

radius [m] radius [m]

Figure 5. Effect of strength of turbulent dispersion force on the radial gas volume fraction profile.

4.5 Turbulent wall function

Prediction of the velocity profiles can be greatly improved by the introduction of new wall function for boiling flow similar to the suggestion of Ramstorfer [28] used subsequently also by Koncar [29]. The idea is that the presence of bubbles growing on the heated wall forces the liquid into a similar flow pattern as that observed in single phase turbulent flow over

a rough wall. For such flows a modified law of the wall applies as described in standard texts on turbulence. For the application to boiling flows we relate here the length scale scharacterizing the hydrodynamic effect of wall roughness to bubble detachment size and nucleation site density as

s rc NdW .

For the constant of proportionality a value of 1 has been found to give good results.

The comparison of models with the single phase (dashed lines) and boiling flow wall functions (solid lines) to the measured velocity profiles in Fig. 6 clearly shows a better agreement of the variant with boiling flow wall function. The simulation results for the other quantities, gas volume fraction, liquid temperature and bubble size, are not affected significantly by the chosen wall function.

■2*2.5 o

DEBORA1 * * * UGAS Exp

■ UGAS

j Uloud --■-i-■-i-'- —'— —'—

0 0.002 0.004 0.006 0.008 0.01 radius [m]

0 0.002 0.004 0.006 0.008 0.01 radius [m]

Figure 6. Effect of turbulent wall function on the velocity profiles (dashed lines phase wall function, solid lines - boiling flow wall function).

single

4.6 Hydrodynamic boundary condition

Related to the wall function is the question of the appropriate hydrody-namic boundary conditions for liquid and gas phases, which apparently has not received due attention in the literature so far. For adiabatic flows it is reasonable to assume that there is no contact between the gas bubbles and the walls. This makes it immediately obvious that appropriate boundary conditions are a no-slip condition for the liquid velocity and a free-slip condition for the gas velocity. For boiling flows, basically two different scenarios may be envisioned depending on whether the bubbles generated at

the wall remain attached to their nucleation site or are sliding along the wall. In the case of sliding bubbles as described e.g. in [30], the boundary condition for the liquid is free slip on the wall area Aw covered by the bubbles and no-slip on the rest. If the bubbles are immobile on the surface, the appropriate condition for the liquid is no-slip on the whole wall area. Since the gas fraction in the two-fluid model describes only those bubbles that have already left the wall, the appropriate condition for the gas phase is free-slip in both cases.

Unfortunately, for the DEBORA tests no information on the occurrence of sliding bubbles is available. Therefore, in addition to the reference calculations which correspond to the assumption that bubble sliding does not occur, we have performed also calculations where the wall shear stress is multiplied by the fraction of wall area covered by bubbles in the application of the wall function boundary condition. Comparison of the velocity profiles from this calculation (solid lines) to the reference case (dashed lines) in Fig. 7 shows that less agreement with the data is obtained with this modification. Thus we may tentatively conclude that bubble sliding likely does not occur in the DEBORA tests. To make a final judgment on the appropriate hydrodynamic boundary conditions in relation to the occurrence of bubble sliding, of course, experiments would be needed where information on this phenomenon is available.

«k„ "I 1 ! .............-VJIiKt, ........— —

--1--- ; DEBORA1 : ^ » « «UGAS Exp *x \ *

uGAS X

Uliquid --'—i—'—i—'—i—'— —'—

«2.0--------'

S 1.5-

0.002 0.004 0.006 0.008 0.01 radius [m]

DEBORA4

» x * UGAS Exp

"-U„„ "

■ U,

LIQUip

-1-'-1-'-1-'-1-'-1

0.002 0.004 0.006 0.008 0.01 radius [m]

Figure 7. Effect of hydrodynamic boundary condition for the liquid phase on the velocity profiles (dashed - no-slip on whole wall area; solid - free-slip on part of wall area covered by growing bubbles).

5 Conclusion

Boiling at a heated wall has been simulated by an Euler/Euler description of two-phase flow combined with a heat flux partitioning model describing the microscopic phenomena at the wall by empirical correlations. Such an approach was previously used and adjusted to boiling experiments with water at a pressure of several MPa.

We have shown that with a careful recalibration of model parameters application to similar tests using R12 at the DEBORA facility is feasible. The need to carefully assess the validity of frequently used correlations for each specific application must be emphasized. It would be highly desirable to develop correlations with an assured broad range of validity or even better well founded physical descriptions of bubble generation and dynamics from which the required closures can be derived.

Beyond the widespread reference model for subcooled flow boiling we have explored several possibilities to improve quantitative agreement with the data. A MUSIG approach has been shown to provide a good prediction of bubble sizes but this does not bring about better predictions of other observables. Bubble forces significantly affect the distribution of void fraction and more detailed modelling is necessary to improve this. The turbulent dispersion force appears least well-founded in this respect. A turbulent wall function for boiling flow has been related to bubble parameters and shown to improve the prediction of gas velocities. It appears likely that the description of bubble induced turbulence in the bulk according to Sato [21] is also too simplified. This would affect several other aspects of the model, namely the turbulent wall function, the turbulent dispersion force and rates for bubble coalescence and fragmentation. Finally, two variants of hydro-dynamic boundary conditions for the liquid have been considered, a no-slip condition corresponding to bubbled growing attached to their nucleation site and a partial slip condition for sliding bubbles. Differences between both are appreciable but to settle this yet rarely addressed issue, more detailed data are needed.

Overall, our results confirm the great potential of the Euler/Euler two-phase flow and heat flux partitioning models for the simulation of subcooled flow boiling in industrial applications while at the same time highlighting the need for specific model improvements in order to achieve highly accurate quantitative predictions.

Acknowledgement This work was funded by the German Federal Ministry of Education and Research under the contract number 02NUK010A.

Received 23 February 2012

References

[1] Ishii M.: Thermo-fluid Dynamic Theory of Two-phase Flow, Eyrolles, Paris 1975.

[2] KüRüL N., Podowski M.Z.: Multidimensional effects in forced convection sub-cooled boiling. In: Proc. 9th Int. Heat Transfer Conf., Jerusalem, Israel, 1990, Vol. 2, Paper 1-Bo-04.

[3] Kurul N., Podowski M.: On the modeling of multidimensional effects in boiling channels. In: ANS Proc. 27th Natl. Heat Transfer Conf., Minneapolis, MN, USA, 1991, 30.

[4] Krepper E. et al.: Modelling of subcooled boiling - concept, validation and application to fuel assembly design. Nucl. Eng. Des. 237(2007), 716-731.

[5] Bartolomej G.G., Chanturiya V.M.: Experimental study of true void fraction when boiling subcooled water in vertical tubes. Thermal Eng. 14(1967), 123-128 (translated from Teploenergetika 14(1967), 80-83).

[6] KrepperE., Rzehak R.: CFD for subcooled flow boiling: .Simulation of DEBORA experiments. Nucl. Eng. Des. 241(2011), 3851-3866.

[7] Garnier J. et al.: Local measurements onflow boiling of refrigerant 12 in a vertical tube. Multiphase Sci. and Technol. 13(2001), 1-111.

[8] Yao W., Morel C.: Prediction of parameters distribution of upward boiling two-phase flow with two-fluid models. In: Proc. 10th Int. Conf. Nuclear Engineering, Arlington, Virginia, USA, 2002, Paper ICONE-10-22463.

[9] Yao W., Morel C.: Volumetric interfacial area prediction in upward bubbly two-phase flow. Int. J. Heat Mass Transfer 47(2004), 307-328.

[10] Boucker M. et al.: Towards the prediction Of local thermal-hydraulics in real PWR core conditions using Neptune_ CFD_ software. Workshop on Modeling and Measurements of Two-Phase Flows and Heat Transfer in Nuclear Fuel Assemblies, KTH Stockholm, Sweden, 2006.

[11] Morel C., Lavieville J.M.: Modeling of multisize bubbly flow and application to the simulation of boiling flows with the NEPTUNE CFD code. Science and Technology of Nuclear Installations (2009), 953527.

[12] Kader B.A.: Temperature and concentration profiles in fully turbulent boundary layers. Int. J. Heat Mass Transfer 24(1981), 1541-1544.

[13] Wintterle T.: Development of a Numerical Boundary Condition for the ¡Simulation of Nucleate Boiling at Heated Walls. PhD. thesis, University of Stuttgart, 2004.

[14] Mikic B.B., Rohsenow W.M.: A new correlation of pool-boiling data including the fact of heating surface characteristics. Trans. ASME J. Heat Transfer 91(1969), 245-250.

[15] Ranz W.E., Marshall W.R.: Evaporation from drops. Chem. Eng. Prog. 48(1952), 141-146.

[16] Anglart H. et al.: CFD prediction of flow and phase distribution in fuel assemblies with spacers. Nucl. Eng. Des. 177(1997), 215-228.

[17] Ishii M., Zuber N.: Drag coefficient and relative velocity in bubbly, droplet or particulate flows. AIChE J. 25(1979), 843-855.

[18] Tomiyama A. et al.: Transverse migration of single bubbles in simple shear flows. Chem. Eng. Sci. 57(2002), 1849-1858.

[19] Burns A.D. et al.: The Favre averaged drag model for turbulence dispersion in Eulerian multi-phase flows. In: Proc. 5th Int. Conf. on Multiphase Flow, Yokohama, Japan, 2004, Paper 392.

[20] Menter F.: Two-equation eddy-viscosity turbulence models for engineering applications. AIAA J. 32(1994), 1598-1605.

[21] Sato Y. et al.: Momentum and heat transfer in two-phase bubble flow-I. Int. J. Multiphase Flow 7(1981), 167-177.

[22] Tolubinsky V.I., Kostanchuk D.M.: Vapour bubbles growth rate and heat transfer intensity at subcooled water boiling. Heat Transfer 1970. In: Proc. 4th Int. Heat Transfer Conf., Paris, France, 1970, Vol. 5, Paper B-2.8.

[23] Kolev N.I.: Uniqueness of the elementary physics driving heterogeneous nucleate boiling and flashing. Nucl. Eng. Technol. 38(2006), 175-184.

[24] Krepper E. et al.: The inhomogeneous MUSIG model for the simulation of poly-dispersed flows. Nucl. Eng. Des. 238(2008), 1690-1702.

[25] Lucas D. et al.: Condensation of steam bubbles injected into subcooled water. In: Proc. 13th Int. Topical Meeting on Nuclear Reactor Thermal Hydraulics (NURETH-13), Kanazawa City, Japan, 2009, Paper N13P1097.

[26] Prince M.J., Blanch H.W.: Bubble coalescence and break-up in air-sparged bubble columns. AIChE J. 36(1990), 1485-1499.

[27] Luo H., Svendsen H.F.: Theoretical model for drop and bubble breakup in turbulent dispersions. AIChE J. 42(1996), 1225-1233.

[28] Ramstorfer F. et al.: Modelling of the near-wall liquid velocity field in subcooled boiling flow. In: Proc. ASME Summer Heat Transfer Conf., San Francisco, California, USA, 2005, Paper HT2005-72182.

[29] Koncar B., Krepper E.: CFD simulation of convective flow boiling of refrigerant in a vertical annulus. Nucl. Eng. Des. 238(2008), 693-706.

[30] Klausner J. et al.: Vapor bubble departure in forced convection boiling. Int. J. Heat Mass Transfer 36(1993), 651-662.

[31] ANSYS CFX-Solver Theory Guide, Release 12.1. ANSYS Inc., 2009.