Group Title: 7th International Conference on Multiphase Flow - ICMF 2010 Proceedings
Title: 17.7.3 - Assessment of the Heat Transfer Model and Turbulent Wall Functions for CFD Simulations of Subcooled and Saturated Boiling
ALL VOLUMES CITATION THUMBNAILS PAGE IMAGE ZOOMABLE
Full Citation
STANDARD VIEW MARC VIEW
Permanent Link: http://ufdc.ufl.edu/UF00102023/00434
 Material Information
Title: 17.7.3 - Assessment of the Heat Transfer Model and Turbulent Wall Functions for CFD Simulations of Subcooled and Saturated Boiling Boiling
Series Title: 7th International Conference on Multiphase Flow - ICMF 2010 Proceedings
Physical Description: Conference Papers
Creator: Prabhudharwadkar, D.
de Bertodano, M.L.
Buchanan Jr., J.R.
Publisher: International Conference on Multiphase Flow (ICMF)
Publication Date: June 4, 2010
 Subjects
Subject: two-fluid model
subcooled boiling
bubbly flow
turbulence wall functions
 Notes
Abstract: This paper describes the details of validation of the two-fluid model for subcooled boiling simulations using the CFD code CFX. Since subcooled boiling occurs at a superheated wall placed in a subcooled liquid, most of the heat and mass transfer takes place close to the wall. Therefore, this study was focused on the assessment of two important aspects of near wall predictions, viz., the wall heat transfer model and the turbulent wall functions. Validation was performed using the state-ofthe- art multidimensional experimental data available in the literature. The wall heat transfer model is based on splitting the wall heat flux into three components, viz., the single phase convection, the evaporation at the bubble surface and the quenching effect at the heated wall after the bubble departure. The evaporation and quenching components use the following three closure parameters: 1. Bubble Nucleation Site Density, 2. Bubble Departure Diameter, and 3. Bubble Departure Frequency. It is difficult to test three separate uncertainties simultaneously and hence only those data sets where the bubble diameter has been measured have been selected. The turbulent wall functions are used to prescribe the near wall boundary conditions for momentum and turbulence quantities and hence influence the near wall velocity and turbulence distribution. The state-of-the-art two-fluid model uses a two-phase k-ε turbulence model (Lopez de Bertodano et al., 1994) with the standard wall function (Launder and Spalding, 1974). Although the law-of-the-wall is universal for single phase flows, the same is not true when the boundary layer contains a two-phase bubbly mixture. Previous experiments (Marie et al., 1997) have shown that the velocity profile still follows a logarithmic profile but the slope and the intercept constant vary with the concentration of dispersed fluid and the relative magnitude of the shear and buoyant forces. A modification of the law of the wall coefficients has been suggested by Marie et al. (1997) based on the experimental data obtained from adiabatic air-water bubbly flow over a vertical flat plate. The wall law coefficients of CFX were modified using the two-phase wall function model and significant improvements were noticed in the turbulent parameter predictions.
General Note: The International Conference on Multiphase Flow (ICMF) first was held in Tsukuba, Japan in 1991 and the second ICMF took place in Kyoto, Japan in 1995. During this conference, it was decided to establish an International Governing Board which oversees the major aspects of the conference and makes decisions about future conference locations. Due to the great importance of the field, it was furthermore decided to hold the conference every three years successively in Asia including Australia, Europe including Africa, Russia and the Near East and America. Hence, ICMF 1998 was held in Lyon, France, ICMF 2001 in New Orleans, USA, ICMF 2004 in Yokohama, Japan, and ICMF 2007 in Leipzig, Germany. ICMF-2010 is devoted to all aspects of Multiphase Flow. Researchers from all over the world gathered in order to introduce their recent advances in the field and thereby promote the exchange of new ideas, results and techniques. The conference is a key event in Multiphase Flow and supports the advancement of science in this very important field. The major research topics relevant for the conference are as follows: Bio-Fluid Dynamics; Boiling; Bubbly Flows; Cavitation; Colloidal and Suspension Dynamics; Collision, Agglomeration and Breakup; Computational Techniques for Multiphase Flows; Droplet Flows; Environmental and Geophysical Flows; Experimental Methods for Multiphase Flows; Fluidized and Circulating Fluidized Beds; Fluid Structure Interactions; Granular Media; Industrial Applications; Instabilities; Interfacial Flows; Micro and Nano-Scale Multiphase Flows; Microgravity in Two-Phase Flow; Multiphase Flows with Heat and Mass Transfer; Non-Newtonian Multiphase Flows; Particle-Laden Flows; Particle, Bubble and Drop Dynamics; Reactive Multiphase Flows
 Record Information
Bibliographic ID: UF00102023
Volume ID: VID00434
Source Institution: University of Florida
Holding Location: University of Florida
Rights Management: All rights reserved by the source institution and holding location.
Resource Identifier: 1773-Prabhudharwadkar-ICMF2010.pdf

Full Text



7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30-June 4, 2010





Assessment of the Heat Transfer Model and Turbulent Wall Functions for CFD
Simulations of Subcooled and Saturated Boiling


Deoras Prabhudharwadkar't, Martin Lopez de Bertodano1'*, John Buchanan Jr.2

'Purdue University, School of Nuclear Engineering, 400 Central Drive, West Lafayette, IN-47906, USA
tdeoras@purdue.edu, *bertodan @purdue.edu

2Bechtel Marine Propulsion Corporation, West Mifflin PA-15122, USA
buchanan@bettis.gov

Keywords: Two-Fluid Model, Subcooled boiling, Bubbly Flow, Turbulence Wall Functions

Abstract

This paper describes the details of validation of the two-fluid model for subcooled boiling simulations using the CFD code
CFX. Since subcooled boiling occurs at a superheated wall placed in a subcooled liquid, most of the heat and mass transfer
takes place close to the wall. Therefore, this study was focused on the assessment of two important aspects of near wall
predictions, viz., the wall heat transfer model and the turbulent wall functions. Validation was performed using the state-of-
the-art multidimensional experimental data available in the literature.
The wall heat transfer model is based on splitting the wall heat flux into three components, viz., the single phase
convection, the evaporation at the bubble surface and the quenching effect at the heated wall after the bubble departure. The
evaporation and quenching components use the following three closure parameters: 1. Bubble Nucleation Site Density, 2.
Bubble Departure Diameter, and 3. Bubble Departure Frequency. It is difficult to test three separate uncertainties
simultaneously and hence only those data sets where the bubble diameter has been measured have been selected.
The turbulent wall functions are used to prescribe the near wall boundary conditions for momentum and turbulence
quantities and hence influence the near wall velocity and turbulence distribution. The state-of-the-art two-fluid model uses a
two-phase k-e turbulence model (Lopez de Bertodano et al., 1994) with the standard wall function (Launder and Spalding,
1974).
Although the law-of-the-wall is universal for single phase flows, the same is not true when the boundary layer contains a
two-phase bubbly mixture. Previous experiments (Marie et al., 1997) have shown that the velocity profile still follows a
logarithmic profile but the slope and the intercept constant vary with the concentration of dispersed fluid and the relative
magnitude of the shear and buoyant forces. A modification of the law of the wall coefficients has been suggested by Marie et
al. (1997) based on the experimental data obtained from adiabatic air-water bubbly flow over a vertical flat plate. The wall
law coefficients of CFX were modified using the two-phase wall function model and significant improvements were noticed
in the turbulent parameter predictions.


Introduction

The prediction of boiling flows depends on accurate
representation of the heat and mass transfer processes at
the wall where the vapor is formed. The accuracy of the
model prediction depends on how well the following
three independent parameters are estimated: 1. Bubble
Nucleation Site Density, 2. Bubble Departure Diameter,
and 3. Bubble Departure Frequency. These parameters
depend on complex physical phenomena that are not
completely understood. Different closure relations exist
in the literature for each of the above parameters which
are discussed in detail in the next section.


One of the goals of the present study is to perform an
independent assessment of each of these parameters from
the perspective of CFD implementation. However, it is
difficult to test three separate uncertainties simultaneously
and hence only those data sets have been selected where the
bubble diameter has been measured. A review of the
departure frequency models was done to select the one best
suited for the present data sets. Then, using the known
departure bubble size and the chosen departure frequency
model, two well known models for the nucleation site
density were implemented and the results were compared.
The heat transfer at the wall also depends on the turbulence
intensity near the wall and this should be correctly
predicted as well. The wall function approach that is









usually used with the Standard k-epsilon model to obtain
the near wall turbulence quantities using coarser meshes
needs a correction that accounts for the effect of presence
of a two-phase mixture. Such an approach has been
implemented successfully for adiabatic flows (Marie et
al., 1997) and this has been extended for boiling flows
here.
The two-fluid model used here was previously tested for
adiabatic flows where the interfacial momentum transfer
terms were validated (Prabhudharwadkar et al., 2009).
The energy equation along with the interfacial heat and
mass transfer correlations were incorporated into the
model. The next section describes in detail the two-fluid
model used in this study.

Nomenclature

ai Interfacial Area Concentration
CD Drag Coefficient
CVM Coefficient of Virtual Mass (0.5)
Db Bubble Diameter
Ddep Bubble departure diameter
ek Phase enthalpy
fdep Bubble departure frequency
g Gravitational Acceleration
hi, Latent heat of evaporation
jL Superficial Velocity of Liquid
jG Superficial Velocity of Gas
k Turbulence Kinetic Energy
Ma Interfacial Momentum Source
Amk Mass Source
nste Nucleation site density
q" Heat Flux
v-k Velocity (vector)
v-R Relative Velocity (vector)
T Temperature
Greek Symbols
ka Volume Fraction of kh phase
a Void Fraction
F Dissipation of Turbulence Kinetic Energy
Pk Density
S Turbulent (Reynolds) Stress
S Dynamic viscosity
A Thermal Conductivity


Superscript
D Drag
TD Turbulent Diffusion
L Lift
W Wall
Subscript
k Phase identifier
1, 1, f, L Continuous (Liquid) Phase
2, g, G Dispersed (Gas) Phase
dep Departure
t Turbulent


7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30-June 4, 2010



Two-fluid Model with Heat and Mass Transfer

For diabetic flows, the ensemble averaged two-fluid model
equations (Ishii and Hibiki, 2006) governing the motion of
each phase has the following form:

Mass Conservation:


ak Pk + V.OkPkV = Amk
d t


Momentum conservation:

dv = akvPk + RY)
,Pk + V ak Pk k k k--VP + k +Z=k
at


+ apk g_ + Mt + ,,, ,v

(2)


Energy conservation:


z, z, Re
d ape + V-a ,U,e,= -V.k c+ qR)


aL H (3)


qkk ktk V, (4)
q(k = -kVTk qk kVek

In the present calculations of subcooled boiling, the vapor
generated was assumed to be at saturated temperature and
hence was treated as isothermal phase.

Interfacial momentum transfer
The interfacial momentum transfer force comprises of force
terms due to drag, turbulent diffusion, lift and a wall force:


M =MD L TD W
Mki Mi +M +M +M
mki -- ki ki ki


The drag force of the bubbles is given by,

D 3 C R R(6)
M C-ap1 D VR (6)
4 D,

where the drag coefficient, CD, is given by the Ishii-Zuber
correlation (1979). Db is the Sauter mean diameter of the
bubble field and the relative velocity is given by
VR =V2 -)1

The wall force model accounts for the effect that keeps the
centers of the bubbles no closer than approximately one
bubble radius from the wall. This force is important when
the lift force is present. It is given by (Antal et al., 1991):









M2 =C ,C (. -. n,

Ca= min 0 c ,, (7)
Db Ywall / ]

where n and ywal are the normal vector and the
distance from the wall respectively. The CFX values for
the coefficients are: Cw = -0.01 and Cw2= 0.05. In a
previous study (Prabhudharwadkar et al., 2009) it was
found that the Antal's CFX Default coefficients work
well for adiabatic bubbly flows through vertical pipe.

The turbulent diffusion force is closed using the Lopez de
Bertodano Model (Lopez de Bertodano, 1991):

TD
M = CTDPk P Va (8)

The turbulent diffusion coefficient used for the present
calculations is CTD = 0.25.


The lift force is given by (Auton (1987)):

M = -CL, pla 2- v)x(vxv,) (9)

The lift coefficient used for the calculations is CL = 0.1
(Lopez de Bertodano, 1991). This value is of the same
order as Tomiyama's experimental value (2002) for a
bubble in Couette flow, CL = 0.288.


Interfacial heat transfer
The interface to liquid heat transfer (k=l) is expressed as,

qj, = h -,a,(Ta,, T,), (10)


where, h1i is the heat transfer coefficient between the
liquid and the interface, closed as follows,

h Nu (11)
d

In the above equation, Al is the liquid thermal
conductivity and d, is a length scale which is assumed
equal to the bubble diameter. The Nusselt number (Nu)
in the above closure equation is given by following
expression (Ranz and Marshall, 1952):

Nu = 2 + 0.6Re5 Pr033 (12)


Reynolds number (Reb) used in the above closure
defined as,

Reb =
vl
The interfacial mass transfer is given as,


7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30-June 4, 2010


h,(,-T,)a h, (Tl -T,)a,
cond,L max 0 eap, = in -, 0

(14)
where, Hig is the latent heat of phase change.


Wall heat transfer
Wall heat transfer is the most important aspect of subcooled
boiling as it provides the sources for energy and phase mass
balance equations. Most wall heat transfer models are based
on partitioning the wall heat flux into three components.
These were implemented in the multidimensional CFD
calculations by Kurul and Podowski (1990). The three
components of the wall heat flux are -
a) Evaporation heat flux (qe) The part of heat flux
utilized in formation of vapor at the wall,
b) Forced convection to the liquid (qc),
c) Wall quenching by liquid phase (transient conduction)
(qq) This accounts for the heat transfer to the
subcooled liquid that replaces the detached bubble at
the wall.

The wall surface is assumed to be split into two parts (A1,
A2) each under the influence of one phase. Fraction A2 is
influenced by the vapour bubbles formed on the wall and
participates in the evaporation and quenching heat transfer.
Fraction A1 is the remaining part of the wall surface,
(AI=1-A2) and participates in the convective heat transfer
to the liquid. A1 and A2 are related to the nucleation site
density per unit wall area (n,i,) and to the influence area of
a single bubble forming at the wall nucleation site. The
Kurul-Podowski model assumes, that the diameter of the
bubble influence zone is twice as big as the bubble
departure diameter (a = 2)


A2 = (a Dde)2 nse A = 1- A2
4


Evaporation Heat Flux:
The evaporation heat flux is obtained from the rate of vapor
generation at the wall which is given as a product of mass
of a bubble at detachment, the detachment (departure)
frequency and the nucleation site density,


SD3
, de PG fdep ns1te

q, = fm A


es is
The specification of evaporation heat flux using above
equation requires closure for bubble departure diameter,
(13) bubble departure frequency and the site density. The default
available model in CFX uses following closure relations:
The bubble departure diameter is given by the empirical
correlation of Tolubinsky and Kostanchuk (1970), the
Nucleation Site Density is obtained using the correlations of
Lemmart and Chawla (1977) and the departure frequency is





7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30-June 4, 2010


closed using the simplest available expression (Cole,
1960), which uses a characteristic bubble velocity
(terminal velocity of bubble rise) and a characteristic
bubble size (departure diameter).
The frequency of bubble departure obtained by Cole
(1960) was derived under the assumption of bubble
moving with terminal velocity once detached from the
surface. This is a fair assumption for pool boiling
scenario at low heat fluxes where bubble detachment is
hydrodynamically governed, but may not be accurate for
flow boiling at high heat fluxes where the
thermodynamics governs the bubble growth and
detachment. The CFX default closure relations are given
below:

Ddep = Dteexp(- (T Tiq )/Tt),
Dr = 0.6 mm, Tef =45K (18)


nsite = n/ ((LTwa Ts ) / ATrs )' 805
n, = 7.9384x105 m2, ATre =10K (19)


= 4gAp ..
d 3DdepPliq

Quenching Heat Flux:
As mentioned previously, the quenching heat flux is the
component of wall heat flux utilized to heat the cold
liquid replacing the detached bubble adjacent to the
heated wall. In order to evaluate this component, Mikic
and Rohsenow (1969) used an analytical approach
starting with transient conduction in semi-infinite
medium with heated wall being the only boundary where
temperature is specified.

DT(x,t) la2T(x,t)
t )=a ; T(x,o)=TF,T(o,t)=T,T(oo,t)=T
dt dx
(21)

Kurul and Podowski (1991) made an assumption that
quenching occurs between detachment of one bubble and
appearance of next bubble (nucleation), and this time
was assumed to be 80% of the detachment period. The
final form of quenching heat transfer closure of Kurul-
Podowski model is as follows:


S2A f t=0.8 (22)
rzra f
q = A2 h(T T ) (23)

Single Phase Convection Heat Flux:
This is evaluated under the standard assumption of a
logarithmic temperature profile across the turbulent
boundary layer (Kader (1981), ANSYS CFX Solver
Theory Guide ,"2 II .,1).


qw,lph = + w w Thql
T+


T+ =Pr, y er +(2.121n(y*)+fp)el/F,

,=f(Pr),F =/f(y),

u =C kY2 y p u* Ay/p,
(25)
where, Ay is the distance of the wall adjacent node from the
wall, ki is the liquid turbulence kinetic energy. In the above
equation fl is a function of liquid Prandtl number and F is a
function of y (ANSYS CFX Solver Theory Guide i'i I
The convective heat flux component through the liquid area
fraction is thus given as,


qc =Alfhc( T,q ), h = -C
T+


It should be noted here that all the wall heat flux model
correlations above use a liquid temperature (Tliq). The
single phase heat transfer correlation above uses the wall
adjacent node temperature. However, in the quenching heat
transfer correlation which was based on a one-dimensional
model, the liquid temperature refers to the bulk mean
temperature. As a good estimate, CFX approximates this
temperature with the temperature at a fixed y (250). This is
done in order to have a mesh size independent evaluation of
wall heat flux partitions.

Review of the Bubble Departure Frequency and
Nucleation Site Density models:
Situ et al. 21"') recently reviewed all the available
departure frequency models and proposed a correlation for
bubble departure frequency based on experimental data of
subcooled boiling flow for wide ranges of pressure and heat
flux (Basu et al. ,'1ll" ), Thorncroft et al. (1998) and Situ
2I'i1n,4 The departure frequency was correlated to the
boiling heat flux as follows:

N = 10.7N .634 (27)

Nfd is a non-dimensional representation of frequency:
Nfd fdep ,
Nfd--
aL,q
where, aLiq is the liquid thermal diffusivity.
NqNB is the dimensionless nucleate boiling heat flux
obtained from Chen's correlation.


qNB Ddep
NqNB -NB D,
aLiqPg hf
q NB hNB (T Tsa,,, )
S(0 79 045 049
hN, = S(0.00122) o1029h24~n


(29)

(30)

T T )024Ap 075

(31)





7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30-June 4, 2010


Ap = p(T)- p(Ts,) (32)

S is a suppression factor which is described in detail in
the next section on site density (equation (41)).

However, when the frequency data of Thorncroft et al.
(1998), Basu et al. (2005) and Situ iI,-1 i), was compared
with the prediction of Situ et al.'s correlation and Cole's
model, the accuracy of the models was found comparable
as shown in the figure below (Fig. 1 (a) and (b)).


(a)
Nfd (experiment)


Nfd (experiment)


Figure 1: Comparison of the bubble departure frequency
models

The data scatter has been contained by the dotted lines in
the figures above which are spaced by an order of
magnitude from the correlation. The data scatter about
the prediction is comparable and hence selection of the
Cole's model (default CFX model) was reasonable for
the present problem.

A correlation for the nucleation site density has been
proposed by Hibiki and Ishii (2003) which is valid over
1-198 bar and for most practical combinations of fluid
and surface material (e.g., Water-Stainless steel, Water-
Copper, Water-Zr-4, R113-Nichrome etc.). This


correlation is essentially an improvement of the previously
established correlation of Kocamustaffaogullari and Ishii
(1983) validated against water data for 0-50 bar. The
Kocamustaffaogullari-Ishii (K-I) correlation relates the site
density to the wall superheat as follows:

n,,e = R44 fR ), (33)
where,
nsie=nte Ddep (34)
R = Rc/(Dde /2) (35)

R, is the critical radius of the surface cavity which
represents a minimum cavity size which can be activated at
a given wall superheat. R, is given as,


2- a(1+ (p, p)) /P
Rc t-.1(, AT /(RTgsa, ))-1)

Under the following conditions:


Pg << p, hg, (T, Ts,) (RTTa,)<<1,


R, can be simplified as,


Rc = 2oTst,, (, hf ATs,,,).


The density dependent parameter in equation (33) is given
as,


f(p*)= 2.157x10 7p*-32(1+0.0049p*)413
P = (-P, )/P,


The effective wall superheat in equation (36) is usually less
than the actual wall superheat. The reason for this is as
follows. A bubble nucleated at the wall grows through a
liquid film adjacent to the wall where a considerably high
temperature gradient exists. So, in reality, it experiences a
lower mean superheat than the wall superheat. In case of
pool boiling, the difference is not significant and the
superheat based on the wall temperature can be taken as the
effective superheat. In case of forced convective boiling,
there exists a steeper temperature gradient in the near wall
field due to thermal boundary layer. Hence, the effective
superheat experienced by the bubble is smaller than the
actual wall superheat and it is given as,


The multiplier S in the above equation is the superheat
suppression factor and is given as (Chen, 1966),


S=1/( +1.5x105 Re1,)


The two-phase Reynolds number is given as,


Re, = [G (1- x)dh /, ]F125


Situ et al.(2008) Correlation Prediction
(Rel


AT,,,, = SAT,, where, AT, = T7 T1
'" sate or"-a at w at





7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30-June 4, 2010


RESULTS


Validation of the Heat and Mass Transfer Models


The Martinelli parameter can be approximated as,


X=,,(= 9 'g 5, f U (44)
Sx A ) p [g)

Thus, equations (33)-(44) represent the complete model
of Kocamustafaogullari-Ishii (1983), which has been
incorporated in CFX-12 for the present study.

Turbulence Transport
The closure for the Reynolds stresses in equation (2) is
based on the two-phase k-e model developed by Lopez de
Bertodano et al. (1994). This model assumes that the
shear induced and bubble induced turbulent stresses are
added together:

Re Re) +(Re) (45)
=k =k SI =k BI

The resulting expression for the diffusivity of
momentum (effective viscosity) of the liquid phase is:
k2 -
Va C -+CDBaDb R

(46)
where, the first term on the RHS corresponds to the k-e
model for the shear induced turbulence viscosity and the
second term corresponds to Sato's model (1981) for the
bubble induced turbulence viscosity. The coefficient C, =
0.09 is the standard value according to the k-e model. A
value of 0.6 is recommended for CDB. It is important to
note that the k-E model transport equations for the liquid
phase are solved together with the continuity and
momentum equations (Eqs. (1) and (2)). The standard
coefficients of the k-s model are left untouched.
Standard Wall Functions (Launder and Spalding (1974))
are used to model turbulence quantities near the wall
which avoids need to resolve the complete boundary
layer. Turbulence kinetic energy dissipation rate
(epsilon) is specified using an algebraic closure
assuming turbulence production due to shear is balanced
by dissipation. Shear stress at the wall is connected to the
velocity in the fully turbulent boundary layer
(computational node next to the wall) using a
logarithmic law of wall which is described in detail later.
This completes the description of the two-fluid model
and all the closure relations used in the present study.
The results of the comparison of two previously
mentioned nucleation site density model are discussed in
the subsequent section. The results reported here are
restricted to the prediction of vapor void fraction
distribution which of primary interest in this study.


1. Subcooled boiling of R12 at 15-27 bar (Morel et al
(2003))
Morel et al. (2003) performed subcooled boiling
experiments in a vertical pipe of 19.2 mm internal diameter
and a length of 5 m. The pipe had a heated section of 3.5 m
preceded by and followed by unheated lengths of 1 m and
0.5 m respectively. The input parameters for the four tests
reported were as summarized in Table 1. The liquid-vapor
density ratio in these experiments corresponds to that for
water-steam at 95-150 bar. Dual sensor fiber optics probes
were used to measure void fraction

Table 1: Simulation input parameters for tests tpl and tp6
Parameter Test1 Test2 Test3 Test4
(Deb5) (Deb6) (Debl3) (DeblO)
Mass Flux 1986 1984.9 2980.9 2027.0
(kg/m2)
Pressure 26.15 26.15 26.17 14.59
(bar)
Inlet 18.12 16.11 18.12 23.24
Subcooling
(C)
Heat flux 73.89 73.89 109.42 76.24
(W/m2) __

A 2" wedge of the pipe is used as the domain. The mesh
had 20 radial and 250 axial uniformly spaced elements.
Fig. 2 shows the domain and cross-section of the mesh
used.


Figure 2: Domain and Radial Mesh for R12 boiling
simulation case

All the data sets include radial profiles of bubble diameter.
Note that with Freon at these test pressures, the surface
tension is of the order of 0.001 N/m which results into very


F=1.0 ,forX,, >10
=2.35(0.213+1/X,,)0736,forX,, <10


(











small bubble sizes (i.e., 0.5 mm). For the CFD
calculations the mean of the radial profile value is used
for the bubble diameter in the bulk, whereas, the value at
the wall is used as the departure bubble size (Fig. 3). It
was found during all the cases that a better match with
experimental data is obtained with the coefficient of the
turbulent diffusion force CTD = 0.5 instead of the value
obtained for adiabatic air water flows, CTD = 0.25. This
change may be attributed to bubble diameters in these
cases which are one order of magnitude smaller than that
in atmospheric air-water systems (Prabhudharwadkar et
al., 2009). Therefore there are more interactions of the
bubbles with smaller turbulent eddies.
Fig. 4 shows the prediction of radial variation of void
fraction obtained using two different site density models.


7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30-June 4, 2010


The Lemmart-Chawla (1977) correlation which is used in
the Kurul-Podowski model was found to under-estimate the
void fraction in all the cases. Kocamustafaogullari-Ishii
model (termed as K-I Model in the figures below) predicted
the data well in all the cases. The prediction however
deteriorated at low pressure data (P = 14.59 bar) and this
needs further investigation.


P = 26.15 bar, q" = 73.89 kW/m2, G = 1986 kg/m2 s,
dTsub = 18.12 C

Experiment (Morel et a] (2003))
-Mean




r ^ ^ ^ urr^ ^ U1ur


2 E-04


0 E+00
0 0002 0004 0006 0008 001
r (m)


P = 26.17 bar, q" = 109.42 kW/m2, G = 2980.9 kg/m2 S
dTsub= 18.12 C

Experiment (Morel et al (2003))
-Mean


4 E-04 f at


2 E-04


0 E+00
0 0002 0004 0006 0008
r (m)


P = 26.15 bar, q" = 73.89 kW/m2, G = 1984.9 kg/m2 s,
dTsub = 16.11 C
Experiment (Morel et al (2003))
-Mean


0 E+00
0 0 002 0 004 0 006 0 008 0 01
r (m)

P= 14 59 bar, q" = 76 24 kW/m2, G = 2027 kg/m2 s,
dTsub= 2324 C
1 5E 03
Experiment(Morel et al (2003))
-Mean

1 OE 03 (m*

**%

5 OF 04 -




0 OE+00
0 0 002 0004 0006 0008 001
r (m)


Figure 3: Bubble diameter approximation for the simulations


1 E-03


8 E-04


6 E-04


4 E-04


2 E-04


1 E-03


8 E-04


6 E-04


~*,

~55







7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30-June 4, 2010


P = 26.15 bar, q" = 73.89 kW/m2, G = 1986 kg/m2 s, P = 2.15 bar, q" 73.89 kW/m2, G =1984.9 kg/m2 s,
dTsub = 18.12 OC, D= 0.452 mm, De,= 0.235 mm dTsub =16.11 C, D = 0.498 mm, Ddp= 0.25 mm


P = 26.17 bar, q" = 109.42 kW/m2, G = 2980 kg/m2 s,
dTsub 17.44 C,D= 0.384 mm, Dd= 0.146 mm



-- --SiteDensiy LemmarChawla(1977)





^^ ^


Figure 4: Void Fraction prediction using the two Nucleation Site Density models


2. Subcooled Boiling of R113 at 2.69 bar in an
Annulus (Roy et al (2002))
Roy et al. (2002) performed subcooled boiling
experiments in a vertical annulus of 15.78 mm internal
diameter and 38.02 mm outer diameter. The total
length of the channel was 3.66 m of which the initial
0.91 m was unheated (adiabatic). The pipe had a heated
section of 2.75 m.
Six experiments were reported (designated as tpl-tp6)
and two of them (having minimum and maximum
mass flux, tpl and tp6) have been selected for
simulation purpose here. Radial Profiles of Void
Fraction were measured with a dual-sensor fiber-optic
probe. The input data for the two representative cases is
stated in Table 3:

Table 3: Simulation input parameters for tests tpl and
tp6_
Parameter tp tp6

Mass Flux (kg/m2) 568 784
Pressure (bar) 2.69 2.69
Inlet Subcooling (C) 42.7 50.2
Inner wall heat flux (W/m2) 95000 116000


Similar to previous problem, a 2 section of the annulus
is used as the domain. The mesh had 25 radial and 266
axial uniformly spaced elements. Fig. 5 shows the
domain and cross-section of the mesh used.






Figure 5: Mesh cross-section for R-113 simulations

These experiments had a characteristic near wall
boundary layer of vapor (boiling layer) and the bulk of
the pipe (more than 50%) had only liquid. Bubble
diameter was reported for the test tp6. As the pressure
was same in all the experiments and the range of other
parameters was also narrow, same values were used for
the other simulated case (tpl). The bubble diameter
profile is shown below in Fig. 6.


P = 14 59 bar, q" = 76 24 kW/m2, G = 2027 8 kg/m2 s,
dTb = 23 24 oC, D=0 916 mm, Ddp=0 37mm

SExperment(Morel et a (2003))
Site Density K (198)
--SiteDensty Le a Chawa(1977)







b.. *..*44W0













------ MipnMe ditmeern
S- ---- S-utr mc dinamccr

2 : --n Ie r
S A,




















assumed to be 1 mm and the departure bubble diameter

0.6 mm arrived at from the "most probable" diameter
values reported in the experiments (Fig. 6).
Simulations were also done with the Sauter mean
diameter values in the bulk and at the wall; however
the predictions with most probable values matched the
data well. As stated earlier, the simulations were
performed with the K-I site density model only since
this data set is in the pressure range for which the K-I
model is already found to be best suited through
previous benchmark problems. The void fraction
prediction for the two cases is shown below in Fig. 7.
The width of the boiling layer was well predicted
with the most probable bubble sizes; however the near
wall void fraction was under-predicted in both the
cases.

Two-phase Wall Function Implementation

The velocity profile in the turbulent boundary layer
next to the wall follows a logarithmic profile as given
below:


U =1log(y+)+c
ic


where, U+=U y+ = uy = (48)
w, Y
U, V p
K is the von Karman constant (= 0.41) and C is a log-
layer constant (= 5.2). U is the velocity component
tangential to the wall and y is the distance normal to
the wall.


7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30-June 4, 2010


0.6
Experiment(Roy etal. 2002)(tpl)
0.5 Db = 1mm, Ddep =0.6 mm
S--Db= 1.7 mm, Ddep 1.2 mm
0.4




0.2


0.1


0.2 0.4 0.6
(r-ri)l/ro-ri)


0.8 1


0.2 0.4 0.6 0.8
(r-ri)/(ro-ri)


Figure 7: Void fraction distribution prediction for tpl
and tp6 with most probable and mean diameter values



The above definitions are in context of single phase flow.
Marie et al. (1997) performed experiments where a flat
plate was placed parallel to bubbly air-water flow. Void
fractions and velocity profiles were measured in the
boundary layer next to the wall for different inlet
concentrations of the air. It was found that the velocity
profile still follows a logarithmic law-of-the-wall but
with a modified slope and log layer constant. When the
streamwise momentum equation for liquid phase was
integrated through the fully developed two-phase
boundary layer, following equation was obtained:


(i-a(vL -UL = L u g(a )BL (49)

where, (ULVL) is the turbulent shear stress (Reynolds

stress), a is the average void fraction in the boundary
layer and aO is the average void fraction in the free











stream. YBL is the boundary layer thickness. The second
term on the right side of the above equation represents
the buoyancy due to presence of dispersed phase.
Neglecting the viscous stresses, the above equation can
be simplified for dilute dispersed flows as,


(ULVL) "L -g(a a )BL = 2 (50)

The two-phase log layer equation proposed by Marie et
al. (1997) uses u, as the velocity scale instead of u, in
equations (47-48). When this velocity scale was used
for the law of the wall it was found that the log layer
profiles for all the inlet void fraction were parallel and
the log layer coefficient C increased with the inlet void
fraction. Fig. 8 below shows the results obtained by
Marie et al (1997) for different inlet void fractions after
using a modified velocity scale.


35

30

25 7-



15
S0 a 0 a = o%
10 A- a a = 2% -
V a a =3.5%
5 0 0 a 6..


o10 10' 12 /1 103

Figure 8: Log layer velocity profiles measured by
Marie et al. (1997) for different air void fractions

The log layer constant for two-phase mixture, C*, was
related to the single phase C with the following
relation.

C*= C+yo+ (1/l -1)-(1/Ic)log1 (51)


y +represents the edge of the laminar sub-layer which
is usually about 11 and

/=U (52)
U
Based on the above model, the log layer coefficients
were modified in the CFD code CFX and turbulence
quantities were compared with single phase and two
phase law of the wall. The above mentioned
experimental benchmark of Roy et al (2002) included
measurements for turbulence kinetic energy and
Reynolds stress. A single phase case was first tried to
validate the single phase wall function included in the
CFD code CFX. Fig. 9 shows turbulence kinetic energy


7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30-June 4, 2010

and Fig. 10 shows the turbulent shear stress profile for a
single phase experiment reported by Roy et al (test spl:
Mass Flux = 568 kg/m2, Pressure = 2.69 bar, Inlet
Temperature = 315.85 K, Wall Heat Flux = 16 kW/m2).
The results were satisfactory confirming adequacy of the
mesh size and model accuracy for the prediction of near
wall turbulence quantities.


0005
Experiment Royetal (2002)
0004 CFXResults

0003


0002

0001


0 02 04 06
(r-ri)/(ro-ri)


08 1


Figure 9: Turbulence kinetic energy prediction for
single phase case


0000o



0000


-0 ooo0005


-0 001


Figure 10: Reynolds stress prediction for single phase
case

The two phase log law coefficients were then
implemented in CFX to compare the effective turbulence
kinetic energy and turbulent stress. The effective
turbulence kinetic energy superposes the shear induced
and bubble induced components as follows:

kL, = kL,SI + kL,BI; kBI = 4 avR (53)

Significant improvement was noticed in the prediction
of the turbulence kinetic energy (Fig. 11) for both the
previously reported cases especially near the boiling
boundary layer edge where a single phase wall function
showed a characteristic dip.
The Reynolds stress (obtained using equations (45-46))
was found to be more accurate with the modified wall
law coefficients (Fig. 12).


Experiment (Royetal (2002))(spl)
CFX Results



4 06 0

S 04 06 08






(r-ri)/(ro-ri)







7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30-June 4, 2010


S 0.2 0.4 0. 0 O 1 0 0.2 0.4 0.6 08 1
(r-ri)/(ro-ri) (r-ri)/(ro-ri)
Figure 11: Turbulence kinetic energy prediction for wall-boiling cases tpl and tp6 of Roy et al. (2002)


0002 0.002

0.001 0.001

0 0
0.2 0.4A 0.6 0.8 0.2 04 0.6 0.8

S 1 -0.001
A A
S p Experiment(Roy etal. (2002))(tp5)
-o002 -0002
Experiment(Roy et al. (2002))(tpl) -0.002
Originalwall function ginWall Function
-003 Two-phae wall function -0.003 --Two-phase wall function

-0004 0.004
(r-ri/liro-ri) .00(r-ri/(ro-ri)

Figure 12: Reynolds Stress prediction for wall-boiling cases tpl and tp6 of Roy et al. (2002)

assemblies with spacer. Nuclear Engineering and
Conclusions Design, Vol. 177, 215-228 (1997)


1. The heat flux splitting model for wall-boiling
flows was validated using data sets where bubble
diameters were known.
2. The nucleation site density model by
Kocamustaffaogullari and Ishii (1983) in
combination with the Bubble Departure Frequency
model of Cole (1960) is found to predict the data
well.
3. Predictions were found to deteriorate at lower
pressures of the R12 data; however the results were
still satisfactory.
4. Use of a modified log law coefficients based on the
approach of Marie et al. (1997) improved
prediction of the turbulence quantities near the
wall.

References

Anglart, H, Nylund, O., Kurul, N., Podowski, M.Z.
CFD prediction of flow and phase distribution in fuel


ANSYS CFX Solver Theory Guide, ANSYS CFX
Release 11.0, (2006)

Antal, S.P., Lahey, R.T. Jr., Flaherty, J.E. Analysis of
phase distribution in fully developed laminar bubbly
two-phase flow. International Journal of Multiphase
Flow, Vol. 17, No. 5, 635-652 (1991)

Auton T. R., Hunt J. C. R. and Prud'homme M. The
force exterted on a body in inviscid unsteady non-
uniform rotational flow. Journal of Fluid Mechanics,
Vol. 197, 241-57 (1987)

Basu, N., Warrier G. R., Dhir, V.K. Wall heat flux
partitioning during subcooled flow boiling: Part II-
Model validation. Journal of Heat Transfer, Vol. 127, 73-
93 (2005)

Chen, J.C. Correlation for boiling heat transfer to
saturated fluids in convective flow. I& EC Proc. Des.
Dev., Vol. 5, 322-329 (1966).


Experiment(Roy etal. (2002))(tpl)
-- OriginalWall Function
- -Two-phase Wall Function







* *^ e"V e- - -.


E
0.004


0002










Cole, R. A Photographic Study of Pool Boiling in the
Region of The Critical Heat Flux. AIChE Journal, Vol.
6, 533- 534 (1960)

Egorov, Y, Menter, E Experimental Implementation of
the RPI Wall Boiling Model In CFX-11, Ansys CFX,
Staudenfeldweg 12, 83624 Otterfing, Germany (2004)

Hibiki, T. & Ishii, M. Active nucleation site density in
boiling systems. International Journal of Heat & Mass
Transfer, Vol. 46, 2587-2601 (2003)

Ishii, M. & Hibiki, T., Thermo-fluid dynamics of two-
phase flow, Springer (2006)

Ishii, M. & Zuber N., Drag coefficient and relative
velocity in bubbly, droplet or particulate flows. AIChE
Journal, Vol. 25, No.5, 843-855 (1979)

Kader, B.A. Temperature and concentration profiles in
fully turbulent boundary layers. International Journal of
Heat and Mass Transfer, Vol. 24, Issue 9, 1541-1544
(1981)

Kocamustaffaogullari G. & Ishii M. Interfacial area and
Nucleation Site Density in boiling systems.
International Journal of Heat and Mass Transfer, Vol.
26, 1377-1397 (1983).

Kurul, N. & Podowski, M. Z., Multidimensionnal
effects in forced convection subcooled boiling,
International. Heat Transfer Conference, Jerusalem,
Vol. 1, 21-26 (1990)

Lahey, R.T., Jr., Lopez de Bertodano M. and Jones
O.C., Jr., Phase distribution in complex geometry
conduits, Nuclear Engineering and Design, Vol. 141,
pp. 177-201, 1993.

Launder B.E., and Spalding D.B., The numerical
computation of turbulent flows, Computer Methods in
Applied Mechanics and Engineering, Vol.3, pp. 269-
289, 1974.

Lemmert, M., Chawla, J.M. Influence of flow velocity
on surface boiling heat transfer coefficient, In: Hahne,
E., Grigull, U. (Eds.), Heat Transfer in Boiling,
Academic Press and Hemisphere, 237-247 (1977)

Lopez de Bertodano, M. Turbulent Bubbly Flow in a
Triangular Duct, Ph. D. Thesis, Rensselaer Polytechnic
Institute, Troy New York, 1991.

Lopez de Bertodano, M. A., Lahey, R. T., Jr., and
Jones, O. C., Jr. Development of a k-e model for bubbly
two-phase flow. Journal of Fluids Engineering, Vol.
116, 128-134 (1994)


7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30-June 4, 2010


Mikic, B.B. & Rohsenow, W.M. A new correlation of
pool-boiling data including the fact of heating surface
characteristics, Journal of Heat Transfer, 91, 245-250,
(1969)

Morel C., Yao W. & Bestion D. Three dimensional
modeling of boiling flow for the NEPTUNE code,
Proceedings of the 10th International Topical Meeting on
Nuclear Reactor Thermal Hydraulics, Korea (2003)

Prabhudharwadkar, D.M., Bailey, C.A., Lopez de
Bertodano, M.A. & Buchanan Jr., J.R. Two-Fluid CFD
model of adiabatic air-water upward bubbly flow through
a vertical pipe with a One-Group Interfacial Area
Transport Equation, Paper FEDSM2009-78306,
Proceedings of the ASME 2009 Fluids Engineering
Division Summer Meeting, Vail, Colorado, USA (2009)

Ranz, W.E. & Marshall, W.R. Evaporation from drops:
Parts I&II. Chem. Eng. Prog., Vol. 48, Issue 3, 141-146
(1952)

Roy, R.P, Kang, S., Zarate, J.A., Laporta, A., Turbulent
subcooled boiling flow: Experiments and Simulations.
Journal of Heat Transfer, Vol. 124, 73-93 (2002)

Sato, Y, Sadatomi, M. & Sekoguchi, K. Momentum and
heat Transfer in Two-phase Bubble Flow. Int. J.
Multiphase Flow, Vol. 7, 179-190 (1981).

Situ, R., Experimental and theoretical investigation of
adiabatic bubbly flow and subcooled boling flow in an
annulus, Ph.D. Thesis, Purdue University, USA (2004).

Situ R., Ishii M., Hibiki, T., Tu, J.Y, Yeoh, G.H., Mori,
M. Bubble departure frequency in forced convective
subcooled boiling flow. International J. Heat and Mass
Transfer, Vol. 51, 6268-6282 (2008)

Thorncroft, GE., Klausner, J.E, Mei, R. An
experimental investigation of bubble growth and
detachment in vertical upflow and downflow boiling.
International Journal of Heat and Mass Transfer, Vol.
41, 3857-3871 (1998)

Tolubinsky, V.I. & Kostanchuk, D.M. Vapour bubbles
growth rate and heat transfer intensity at subcooled water
boiling. Proceedings of the 4th International Heat
Transfer Conference, Vol. 5, Paris (Paper No. B-2.8),
(1970)

Tomiyama, A., Tamai, H., Zun, I., Hosokawa, S.
Transverse migration of single bubbles in simple shear
flows. Chemical Engineering Science, Vol. 57, 1849-
1858 (2002)




University of Florida Home Page
© 2004 - 2010 University of Florida George A. Smathers Libraries.
All rights reserved.

Acceptable Use, Copyright, and Disclaimer Statement
Last updated October 10, 2010 - Version 2.9.7 - mvs