Matlab_多项式运算与代数方程求解器_百度文库
两大类热门资源免费畅读
续费一年阅读会员，立省24元！
文档贡献者贡献于
评价文档：
38页免费17页免费46页免费34页免费45页2下载券 9页免费6页免费23页1下载券17页免费2页免费
喜欢此文档的还喜欢15页免费42页免费34页免费12页免费20页免费
Matlab_多项式运算与代数方程求解器|计​算​机​语​言​m​a​t​l​a​b
把文档贴到Blog、BBS或个人站等：
普通尺寸(450*500pix)
较大尺寸(630*500pix)
大小：304.00KB
登录百度文库，专享文档复制特权，财富值每天免费拿！
你可能喜欢第3 章 matlab的符号运算第三次课_百度文库
两大类热门资源免费畅读
续费一年阅读会员，立省24元！
文档贡献者贡献于
评价文档：
21页1下载券11页3下载券21页3下载券46页1下载券53页免费 15页4下载券29页2下载券30页2下载券44页免费28页免费
喜欢此文档的还喜欢41页1下载券21页1下载券41页1下载券60页免费46页1下载券
第3 章 matlab的符号运算第三次课|最​好​的​M​a​t​l​a​b​教​程​,​大​学​课​程​教​学​专​用​。
把文档贴到Blog、BBS或个人站等：
普通尺寸(450*500pix)
较大尺寸(630*500pix)
大小：59.50KB
登录百度文库，专享文档复制特权，财富值每天免费拿！
你可能喜欢Services on DemandArticleIndicatorsCited by SciELO Related linksSharePrint version ISSN J. Braz. Soc. Mech. Sci. & Eng. vol.26 no.3 Rio de Janeiro July/Sept. 2004 http://dx.doi.org/10.-00007
TECHNICAL PAPERS
The phase distribution of upward
co-current bubbly flows in a vertical square channel
A. de M E. S. R F. A.
Multlab -Dept. of Energy, FEM -
Campinas, SP. Brazil,
In this work one shows experimental
data and numerical results of the void fraction distribution in vertical upward
air-water bubbly flows in a square cross-section channel. To measure the void
fraction distribution one used a single wire conductive probe. The averaged
void fraction ranged from 3.3% to 15%; the liquid and the gas superficial velocities
varied from 0.9 m/s to 3.0 m/s and 0.04 m/s to 0.5 m/s, respectively. The experimental
results for the void fraction distribution were compared with numerical calculation
performed by an Eulerian-Eulerian implementation of the Two-Fluid Model. In
this work one performs the turbulence modeling with three approaches: using
an algebraic model, the k-e two-phase model and the
k-e two-phase two-layer model. Comparisons between
the experimental and numerical data revealed, in general, good agreement.
Keywords: Vertical upward
bubbly flow, phase distribution in bubbly flows, bubbly flow in square channel
Introduction
It is usual that two-phase bubbly
flows occur in channels of non-circular cross section. Equipment having these
channels are evaporators used in chemical industries, nuclear light-water reactors,
gas-liquid separators used in petroleum plants, to mention a few. The phase
distribution on the cross-section of channels influences the transport processes.
This explains the vast amount of experimental data, modelling and analysis effort
focusing the phase distribution in vertical bubbly flows in round pipes. Contrasting,
there are not so many reports presenting measurements and numerical calculation
for bubbly flows in non-circular channels. One of the reasons is, probably,
due to the experimental difficulties in measuring the void distribution in sections
that do not have axial symmetry. Other is the fact that bubbly flows in oblong
pipes takes the modelling an its numerical implementation to the limit. The
constitutive equations expressing interfacial transfer are not fully established
and there is not enough experimental data to support turbulent models that are
physically robust. Moreover, the implementation of a multidimensional Two-Fluid
Model requires computational effort, power and time, to be solved with an adequate
level of convergence and resolution. Hence, there are just a few papers presenting
numerical results of bubbly flows in non-circular channels and, moreover, comparing
them with experimental data.
Most of the work done on gas-liquid
flows in non-circular channels focused on the mechanistic modelling of pattern
transitions or approach them with one-dimensional models to calculate pressure
drop or mean void fraction. These were the subjects of Mishima et al. (1994),
Troniewsky and Ulbrich (1984), Hibiki et al. (1994), Mishima et al. (1993) and
Keska and Fernando (1994).
Among the papers presenting measurements
of the void fraction distribution in channels having corners one may mention
Sadatomi et al. (1982), Shiralkar and Lahey (1972) e Moujaes and Dougall (1987).
The former generated a vast data bank on vertical air-water flows in rectangular,
isosceles-triangular and annular channels, directed towards the calculation
of the pressure gradient. The authors mentioned that the void distribution is
one of the most important aspects of two-phase flows in non-circular channels.
Their measurements included the cross-sectional distribution of void fraction
in bubbly flows, which depicted wall and corner peaking. The later depicted
data on void fraction distribution, gas and liquid velocity distribution and
bubble size in narrow rectangular channels using both an optical fiber probe
and a hot film anemometer. Some of their results for the void fraction, however,
were presented in terms of the averaged values over the short and long dimensions
of the channels, respectively. This procedure filtered the local distribution
data and revealed only a limited aspect of the void distribution across the
channel section. Non-filtered void distribution data was depicted only in the
longer section center-plane. Similarly to the work of Shiralkar and Lahey (1972),
who focused on bubbly flows in eccentric annulus, various authors dealt with
more complicated flow geometry, like those taking place in nuclear reactor fuel
rod bundles. For the sake of shortness, they will not be mentioned herein.
Lahey et al. (1983) and Lopes de
Bertodano et al. [1994a, 1994b) measured and modelled air-water vertical bubbly
flows in triangular channels. Using hot film anemometers, single and X-probes,
they measured the void fraction and the instantaneous velocity components all
over the channel cross-section, and took fine grid measurements on the corners.
From these they were able to calculate the mean velocity and the Reynolds stress
components, revealing the turbulent structure of the flow. Again, the void fraction
distribution depicted wall and corner peaking, with values somewhat higher in
the corners.
Hoping to improve the understanding
of the phase distribution in non-circular square channels, this work presents,
discusses and compares experimental and numerical data on vertical upward bubbly
flows taking place in a square cross-section channel. The side of the channel
was 0.034 m. The flow was at near atmospheric pressure (0.3 bar to 0.55 bar)
and ambient temperature (23?C). To measure the void fraction distribution
one used a single wire conductive probe. The averaged void fraction ranged from
3.3% to 15%, for liquid and gas superficial velocities varying from 0.9 m/s
to 3.0 m/s and 0.04 m/s to 0.5 m/s, respectively. The measurements depicted
the bubbles concentrating close to the walls and a strong void concentration
in the channel corners.
To model the flow one used the Two-Fluid
Model. A steady-state isothermal flow with no phase change has been considered.
The phase distribution over the channel section came out from a balance between
the radial pressure force, the lateral lift, the wall force and the turbulence-induced
lateral pressure field. One expects the turbulence playing an important role
on the phase distribution over the channel cross section. For this reason, three
different turbulence models have been tested: an algebraic model, the standard
k-e added with a bubble induced-turbulence term and
the k-e two-layer model. The numerical results and
its comparison with the experimental data are presented and commented. For convenience,
Phoenics, a CFD software that utilizes the finite volume approach, solved the
system of equations.
Nomenclatures
Cpi, CD, Cvm,
Cli = pressure, drag, virtual mass and lateral lift coefficients,
dimensionless
Ce2 = constants in turbulent model
Cm = constants in turbulent model
Cpipa, Cla
= constants in lateral lift equation
D = bubble diameter, (m)
E = constant
g = gravitational field (m/s2)
k = Von Karman's constant (dimensionless)
le = scale length, (m)
MkD = generalized
interfacial force, (N/m3)
= generalized interfacial force, (N/m3)
Ng = phase density function,
dimensionless
p = pressure, (N/m2)
P = kinetic energy production, (J)
r = radial displacement
Re = Reynolds number, dimensionless
t = Time, (s)
u = axial velocity, (m/s)
V = velocity field, vectorial quantity,
We = Weber number, dimensionless
x, y = distances from wall, (m)
x,y,z, = cartesian coordinates,
= local void fraction, dimensionless
= delta de Kroenecker, dimensionless
= dissipation rate, (m2/s3)
= turbulent kinetic energy, (m2/s2)
= density, (kg/m3)
nT = kinematic viscosity, turbulent two-phase
kinematic viscosity, (m2/s)
nBI = shear induced, bubble induced kinematic
viscosity, (m2/s)
= surface tension, (N/m)
sk, se = constants in turbulence model
= bubble radius, (m)
= stress tensor, (N/m2)
= turbulent stress tensor, (N/m2)
= azimuthal position
x = rotacional operator
= gradient operator
& & = indicates mean value
Superscripts and Subscripts
k=G , k=L, ki , w = subscripts indicating
the phases (gas and liquid)
and the property at the interface
and at the wall
+ = superscript indicating
a non-dimensional quantity
D = indicates generalized force
d = indicates drag force
tp = subscript indicating a two-phase
The Experimental Apparatus and
the Measurement Technique
a sketch of the experimental set-up. The air-water bubbly flow took place in
a vertical aluminum channel of square cross-section, whose side was 34.1 mm.
In the mixing chamber the air entered through an inner porous media of cylindrical
shape, positioned along the channel axis. This arrangement produced a fairly
uniform bubbly flow just upstream the chamber in the full range of air and water
flow rates applied. After the mixing chamber the bubbly flow developed along
73 equivalent diameters before reaching the test section. The test section was
made of Plexiglas, exactly matching the cross section of the aluminum channel,
allowing flow visualization.
Ordinary tap water and air were the
working fluids, flowing at nearly atmospheric pressure, from 0,3 barg to 0,55
barg, and ambient temperature of 23?C. A centrifugal pump circulated
the water. The air came from a compressed air tank maintained at constant pressure.
To measure the air and the water flow rates within 2% of overall accuracy, one
used a laminar element and a vortex meter, respectively. The range of water
and air superficial velocities - referred to the channel hydraulic diameter
at the test section conditions - applied in the experiments were jL
= 0.9 m/s to 3.0 m/s and jG = 0.04 m/s to 0.5 m/s, respectively.
The mean void fraction in the test section varied from &a&
= 3.3% to 15%. The bubble size was in the range of 2.5 mm to 3.5 mm, measured
from digital images generated by a fast (1000 frames/sec) movie camera. As the
channel had square cross-section, the mean bubble diameter was taken, sampling
a number of bubbles from various frozen images of the channel lateral view,
by comparison with the channel side.
A fully description of the conductive
probe used to measure the void fraction, including the electronics and the data
acquisition can be found in Dias et al. (2000) and Matos et al. (1999).
An inner insulated copper wire of
120 mm diameter, was the cathode in the electrical
circuit. The anode was the probe stem, i.e., the stainless steel needle that
sustained the inner wire. The wire was cut in such a way that exposed a cross
section area having, roughly, 1,1 x 10-8 m2. When immersed
in water, the exposed tip of the wire established electrical contact with the
probe stem through the water. When the tip was inside a gas bubble the electrical
circuit was open. This is the typical on-off signal generated by a classical
conductive probe, from which one calculates the local void fraction. The void
fraction uncertainty does not follow straightforward from experiments like this,
as non-intrusive direct measurements were not performed for comparisons. However,
it can be acessed from global measurements of liquid and gas superficial velocities,
see Dias et al. for details.
scheme of the assembly used to measure the local void fraction over a curved
surface enclosed by the channel walls. To perform the measurements the probe
was inserted into the channel through a small hole and was driven by a mechanism
that allowed radial, r, and rotational positioning. The probe was then positioned
at 105 points (measuring stations) over this surface. These measuring stations
were then projected onto a channel cross-sectional plane – the measuring reduced
plane – where they could be specified by the coordinates (x, y) at the axial
location z. This resulted in a detailed void fraction mapping, required for
a proper averaging process in a strongly non-uniform distribution. A micrometer
with a digital read-out measured the radial displacement, r, with an accuracy
of & 0,02 mm. The azimuthal position, q, could
be set within & 1?. The fact that the actual measurements were
not performed over a cross-sectional plane is of minor significance since the
flow is developed and the axial distance between the most apart measuring stations
is quite short.
The hardware for the data acquisition
consisted of National Instruments& AT-MIO16, a 12 bits ISA bus
board and a 486 PC computer. The hardware sampled, acquired and stored the signals
generated by the probe at every measuring station, at a frequency of 50 kHz,
during the period Dt of 10 seconds. To calculate
the local void fraction from this raw signal, one applied a threshold level
and converted it to a Dirac-delta function. This function - in other words,
the phase density function Ng(x,y,z,t) - indicates the occurrence
of gas (a bubble) at the measuring station on time t. Thus, the local void fraction,
a(x,y,z), comes from the time averaging of Ng(x,y,z,t),
where Dt is the sampling
period. Details of the data processing technique can be found in Dias et al
(2000). It is worthwhile to mention that using the above-described procedure
the number of digitally acquired points used to calculate the void fraction
distribution over the channel section amounted to 4.2 x 107.
The Mathematical Modelling
One used a rather complete formulation
of the Two-Fluid Model, Ishii (1975) and Drew and Lahey (1979), for the mathematical
representation of gas-liquid bubbly flows inside a vertical channel of square
cross-section. A balance between the radial pressure force, the lateral lift,
the wall force and the turbulence-induced lateral pressure field governs the
phase distribution on the channel section. Interfacial transport in two-phase
flows is embodied in the mass, momentum and the interfacial transfer equations.
In a steady-state isothermal flow with no phase change, these equations reduce,
after some definitions and operations, to Eqs. (2), (3) and (4):
The interfacial momentum transfer
equation links the phases:
The subscript k indicates the phase
(k=G or k=L, gas or liquid) and the subscript ki defines proper values
of the variables at the interface, like the velocity and the pressure of phase
k at the interface. The local void fraction is a,
r is the density, V is the velocity, p is the pressure,
g is the gravitational field, tkRe
is the turbulent stress tensor, MkD is the generalized
interfacial force and Mks stands
for the force due to the surface tension. The viscous tensor in both phases
was considered to vanish when compared with tkRe.
In order to solve this system of
equations for gas-liquid bubbly flows, one has to constitute it with proper
equations for MkD, Mks,
the phasic pressure difference, (pki-pk), tkRe
and other variables. Besides some questioning, the first three terms are reasonably
established. The mathematical representation of the two-phase Reynolds tensor,
tRe, however, is far from being settled.
It plays an important role in shaping the multidimensional structure of bubbly
flows since there is a coupling between the two-phase flow turbulence and the
phase distribution. In a developed flow, for example, the radial force balance
accounts for the radial components of forces embedded in the generalized interfacial
force, the radial pressure field and the turbulence-induced lateral pressure
field. As mentioned by Lee (1988), if the
subtle changes in the void distribution are to be modeled, the non-uniform distribution
and the anisotropy of liquid phase turbulence have to be taken into account.
The Reynolds Stress Tensor
As the gas is much less dense than
the liquid, as is in bubbly flows herein reported, it is enough to constitute
the Reynolds, or turbulent, stress tensor in the liquid phase, tLRe.
Every known model available to constitute it has limitations and drawbacks.
Some are not robust enough to represent the physics of the phenomenon, others
bring additional complexities to the system of equations, failing the convergence
of the numerical solution. If one adds the fact that some disagreement persists
in regard to experimental data provided by different authors, it is clear that
the subject of turbulence modelling in two-phase bubbly flows in far from being
settled. Being aware of this, in this work one accomplished it with three approaches:
using an algebraic model, the k-e two-phase model
and the k-e two-phase two-layer model. So, the Reynolds
stress tensor is written according the Boussinesq eddy-viscosity approximation:
where k is the turbulent kinectic
energy and nT is the two-phase turbulent kinematic viscosity, respectively.
Following Bertodano et all (1994c),
one assumes that the two-phase turbulent viscosity is a superposition of the
viscosity due to the Reynolds stress, the shear-induced viscosity, nSI,
plus a bubble induced viscosity, nBI:
The bubble-induced turbulent viscosity
stems from the mixing length concept,
is equal to 0.69 and d is the mean bubble diameter.
In order to solve the phasic mass
and momentum equations one needs to define closure laws for the shear induced
viscosity. Three approaches were used to calculate nT:
an algebraic relation, mentioned herein as the l-vel model, the standard k-e
model and the k-e two-layer model. The l-vel model
is said not being robust, i.e., it does not represent the physics of the phenomenon.
Specifically, the signal change of the Reynolds stress over the pipe diameter,
originated by the non-centered peak of the velocity profile that experimental
data have shown. It is, however, easy to implement and does not impose complexities,
like additional differential equations that could lead to systems that do not
converge. The k-e is, nowadays, a turbulence reference
model, with a widespread use in single-phase flows. However, the standard k-e
model is referred to not giving a good representation of the Reynolds stress
tensor close to channel walls, i.e, when the Reynolds number is low, besides
requiring a very high numerical resolution due to the steep gradient of dissipation
rate of turbulent energy, e, in these regions.
In order to improve computational efficiency and convergence rates, as well
as to introduce established length-scale distribution in the region near the
walls, one alternative is to use the standard k-e
model in the fully turbulent region and an algebraic expression for e in
the damping layer. This is known as the k-e
two-layer model, the third approach one used to calculate nTtp.
The algebraic model (l-vel).
The algebraic model is set by the
non-dimensional relationship between the distance from the wall, y+,
and velocity, V+. The classical Spalding's (1994) law of the wall
where the constants k (von
Karman) and E are equal to 0.417 and 8.6, respectively, and the non-dimensional
variables are
The derivative of y+ in regard to
V+ is the non-dimensional kinematic viscosity, nSI+,
which, by its turn, is .
Hence, differentiating Eq. (6) in regard to V+ gives rise to the
turbulent viscosity:
Now the shear induced viscosity is
set by the k-e model,
where k is the turbulent kinetic energy, e is
the dissipation rate and the constant Cm is equal
to 0.09. The turbulent kinetic energy and the dissipation rate come from the
simultaneous solution of the two-phase version of the the standard k-e model
[19]:
In these equations, sk
= 1, se = 1.314, sa
= 1, Ce1 = 1.44, Ce2
= 1.92. P is the kinetic energy production term:
The single phase "law of the
wall" equation bridges the fully turbulent flows with the flow near the
wall. The values of k and e near the wall come from:
and y are the wall shear stress and the normal distance from the wall, respectively.
two-layer model.
The strategy when using the k-e
two-layer model, according to Chen (1995), is to distinguish flow regions. The
standard k-e model is applied in the fully turbulent
region. In the layer near the wall, an one-equation model, the k-e
equation combined with an algebraic expression for e
, is used:
In the damping layer near the wall
the eddy viscosity is
The Interfacial Pressure Term
In the liquid phase, the pressure
at the interface, pLi, is
where pL is the bulk
liquid pressure.
This expression results from
the solution of a non-viscous flow around a sphere. In this case the coefficient
Cpi is 0.25. Lance and Bataille (1991) suggested values between 0.5~0.7
for actual bubbly flows. Lahey (1990) used Cpi = 1 in distorted bubbly
flows. The value used in the present work was a function of void fraction, as
suggested by Spalding (1994),
where Cpipa is 0.5. Following
Stuhmiller (1977) one considered that the pressure difference is equal for both
phases. Using the Laplace equation the relationship between the gas and the
liquid pressure at the interface could be stated:
is the surface tension and x is the bubble radius.
The Interfacial Forces
The generalized interfacial
force acting on the phase k, MkD, is usually divided in
the drag force, Mkd, and the non-drag force, Mknd.
The non-drag force accounts for the lateral forces and the virtual mass force.
The drag force arises due to the relative displacement of the bubble. The drag
force on a bubble of diameter D is:
According Lahey (1990), the drag
coefficient, CD, depends on the Reynolds number based on the liquid
properties. For Antal et al. (1991), conversely, the Reynolds number must be
calculated for the gas-liquid mixture properties. In this work one follows the
proposition of Kuo and Wallis (1998). According to these authors, the drag coefficient
where Re is the Reynolds number
and We is the Weber number, both expressed as a function of the bubble diameter
and the bubble relative velocity, (VG-VL). The expressions
for CD are valid for bubbly flows in an infinitum medium. Mari&
(1987) showed that next to the pipe wall the drag force alters. In this work,
however, one does not account for the wall effect on the drag force.
A lateral force, MLL=-MLG,
arises when the bubble displaces through a liquid with a non-uniform velocity
distribution. In an upward bubbly flow it acts to displace the bubble toward
According to Drew and Lahey
(1987), the lift coefficient CL varies between 0.05 and 0.1 in air-water
bubbly flows. More recently their co-authors have suggested greater values:
Antal et al. (1991) used values in the range 0.01 ~ 0.5 and Wang et al. (1987)
adopted 0.5. The correlation used in this work was due to Spalding (1994)]:
One varied CLa in
order to find the most appropriate value fitting the experimental data and found
CLa = 0.5. The resulting values for CL were between 0.22
In bubbly gas-liquid flows the bubble
expansion eventually changes the relative displacement between the gas and the
liquid, giving rise to a virtual mass force. The virtual mass force is, Drew
and Lahey (1987):
As the virtual mass coefficient, Cvm,
relates with displaced volumes, one expected it to depend upon the void fraction,
a. However, various authors have used constant values for Cvm. Drew
and Lahey (1987) suggested 0.5; Kuo and Wallis (1988) used values within the
range 2.0 ~ 3.0, and Lance and Bataille (1991) adopted values ranging from 1.2
Similarly with has been done with
the lift coefficient, the virtual mass coefficient has been adjusted to fit
numerical and experimental data. The values that have assured the best agreement,
i.e., the "ideal values" of CL and Cvm, will
be subject of prospective comments.
The Numerical Method
The computations were performed
in a workstation SUN Ultra10, running the CFD PHOENICS 3.2, the package developed
by Spalding (1994). PHOENICS is a transport equation solver with two-fluid model
capability. It is based on the SIMPLE algorithm. The IPSA algorithm simulated
the interfacial momentum exchange.
A computational grid was set up
with 20x20x50 nodes, corresponding to the coordinates (x,y,z), respectively.
A marching procedure, exploring the parabolic nature of the flow, was used to
solve the conservation equations. The initial condition, at the mixer position,
matched up to uniform velocity field and void f the pressure
was fixed, indirectly, through the setting of the gas density. The calculation
then proceeded along 73 equivalent diameters until reaching the measuring reduced
When working with the k-e turbulence
model, careful attention was given when setting the first control volume for
numerical integration, in order to follow the constraints imposed by the "law
of the wall". With the Reynolds number ranging from 30000 to 40000, the
outermost face of the first control volume was located at 2.5 mm from the wall.
This makes the position of the first measuring station different from the first
numerical volume. Hence, the numerical and the experimental results at the measuring
station closest to the wall, in every xz plane, should not, ideally, be compared.
When working with the algebraic model and the k-e
two-layer model this limitation does not exist: the numerical grid could be
refined to rather small distances from the wall.
Regarding the numerical convergence,
one adjusted the linear under-relaxation parameter, LINRLX, used in PHOENICS.
This parameter varies between 0 and 1 and updates the current variable by a
percentage of the previous value. If LINRLX is equal to one there is no under-
if LINRLX is equal to zero, the variables are not updated between two steps
of the integration process. Due to the strong coupling between the phases, the
set of equations had to be under-relaxed at different ratios. The pressure used
a 0,8 factor, the void fraction was under-relaxed at 0.4; the velocities in
the principal and secondary directions were under-relaxed at 0.7 and 0.01, respectively.
The other flow variables, such as the turbulent kinematic viscosity and the
relative velocity were under relaxed at 0.15. These values might not be the
optimum under-relaxation factor, but assured a steady and converging solution.
To reach this converged solution, every operational condition, established by
the gas and liquid superficial velocities, required 5 hours of computation.
Experimental and Numerical Results
Initial tests were conducted to
access the verticality of the channel. The void fraction distribution over the
channel cross-section, appearing in the contour map in , showed that reasonable vertical positioning had been achieved. The software
Surfer& produced the contour map.
Even using a less refined measuring
grid, these first measurements pointed out to results that occurred in every
operational condition: bubbles concentrating near the wall, specially on the
channel corners. The void peaking near the wall, typical in vertical upward
bubbly flows in round pipes, also occurred in the square cross-section channel.
Definite and more detailed measurements
were performed on half of the test section. Exploring the axis-symmetrical nature
of the flow, the local void fraction was measured in 105 stations for various
operational conditions. The measuring grid was particularly refined in the channel
region close to the wall. Corresponding numerical results for the void fraction
distribution have been computed. The numerical results refer to: (i) the three
turbulent models described in a previous section, and (ii) specific formulations
for the lift and virtual mass forces with "ideal values" for CL
and CVM, as discussed.
a set contour maps showing the experimental and numerical distribution of the
void fraction over the channel cross-section. The maps have been plotted using
a commercial package, Surfer&. Due to the fact that one does
not have control over the interpolation procedure embedded in the package, it
is not possible to state that they are a true representation of the data. However,
they are a good overview of the void fraction distribution over the channel
section and permit a quick, besides rough, comparison between the experimental
and numerical data.
To the left, the outermost map represents
the experimental data. The three subsequent maps to the right stand for the
numerical results, each one representing a turbulence model, k-e two-layer,
standard k-e and algebraic l-vel. Every sub-set stands for a pair of superficial
velocities: the upper, #1, (JL= 2,12 m/s; JG = 0,13 m/s);
in the middle, #2, (JL = 1,33 m/s; JG = 0,14 m/s) and
the lowest, #3, (JL= 0,84 m/s; JG = 0,03 m/s). The respective
mean void fractions are: #1: &a& = 5.8% ; #2:
&a& = 9.5% and #3: &a&
= 3.4%. The numerical results were obtained, after extensive calculations, for
"ideal values" of the lift and the virtual mass coefficients, i.e.,
the ones that fitted the best the void fraction distribution over the channel
Looking at the contour maps one
can state that, in a broad sense, the present implementation of the Two-Fluid
Model was able to represent the phase distribution of an ascending vertical
bubbly flow in a square cross-section channel. Moreover, the k-e two-layer turbulence
model, among those used as constitutive equations, gave the best representation
of the phase distribution in the range of the applied superficial velocities.
The "ideal value" for the lateral lift coefficient depended on the
local void fraction and was in the range CL = 0.22 ~0.48; for the
virtual mass coefficient, was CVM = 2.
The analysis of the Cartesian plots
appearing in ,
provided "de-facto" comparisons between experimental and numerical
data. Every plot in the figures refers to the void fraction distribution in
a transversal y-z plane defined by a x position, each figure relating to a pair
of (JL;JG). The lowest the x value, the closest to the
lateral wall the y-z plane is. The highest x value, x= 17.35 mm, defines the
mid-transversal plane.
When one looks at the void distribution
on the two planes close to the wall, x= 0.8 mm and x = 1.8 mm, it is clear that
the numerical distribution, in general, did not have a good agreement with the
experimental data. In these planes, ning from corner to corner, the experimental
void fractions were, in general, greater than the numerical results. The difference
amplified as the mixture velocity decreased, from
to . Also, the experimental distribution showed subtle
changes the numerical calculations were not able to pair. The statement holds
true no matter the turbulence relation embedded in the Two-Fluid Model and the
magnitude of the gas and liquid superficial velocities. There are several factors,
different in nature, which could explain the discrepancies between experimental
and numerical results.
The experimental distribution lacked
symmetry in some cases, pointing out to experimental difficulties in measuring
the void fraction distribution at locations very close to the lateral wall or
in the channel corners. The probe positioning could be argued as one of these
difficulties. Furthermore, at these locations, where the velocities are lower
and secondary flows exist, the probe might interfere with the flow field. If
one looks at the similitude between the experimentation and the mathematical
representation of the flow, new problems arise. The actual distribution of the
void fraction at locations close to the wall depends on a "physical"
characteristic of the flow, the bubble size distribution function. The Two-Fluid
Model does not have the bubble size as a constraint when calculating the void
fraction distribution. In other words, in real world the gas volume has a finite
dimension – the probe only "sees" bubbles greater than a certain size,
for example, while the Two-Fluid Model deals with the mixture as a continuum.
Hence, it might not make sense to compare experimental and numerical results
in locations close to the wall if the calculation grid was smaller than the
bubble size, for example. This was the case of the measuring stations on the
first transversal planes, which had x = 0.8 mm and x = 1.8 mm.
When one analyses the experimental
and numerical void distribution on inner planes, from x= 3.5 mm to x = 17.35
mm, a different figure appeared. Good agreement was achieved, no matter the
magnitude of the mixture velocity or the mean void fraction. Among the turbulence
models used, the k-e two-layer gave, consistently,
the best results. On inner planes the void fraction profile was quite flat in
the channel center region and the numerical results delivered by the Two-Fluid
Model with the k-e two-layer model embedded matched
up to it. Getting closer to the walls the experimental and numerical profiles
were somewhat detached. However, the numerical results delivered by the Two-Fluid
Model with the k-e two-layer model embedded had the
highest gradient. Thus, if not agreed with the experimental data in shorter
distances from the wall, at least had a similar trend.
The numerical results of the Two-Fluid
Model with the standard k-e and the algebraic l-vel
models embedded did not pair the experimental profiles as close as the solution
provided by the k-e two-layer model. The numerical
results were higher than the experimental ones in the channel center region,
and were lower close to the wall.
Conclusion
A vast number of measurements of
the void fraction distribution in an air-water vertical bubbly flow in a square
cross section channel has been made and disclosed. To perform the measurements
a single wire conductive probe has been used. The experimental data revealed
that the void fraction profiles presents the wall peaking that is typical in
upward bubbly flows in round pipes. Moreover, the void concentration in the
pipe corners was particularly high. The void fraction profile has been calculated
using a rather complete implementation of the Two-Fluid Model, constituted with
three turbulence models: an algebraic l-vel model, the standard k-e and the
k-e two-layer models.
The Two-Fluid Model with the k-e two-layer
model embedded gave the best representation of the void fraction distribution.
Close to the walls and on the channel corners the experimental and numerical
profiles were somewhat detached, in some cases. However, the solution provided
by the k-e two-layer model, having the highest void
gradient next to the wall, showed the same trend as the experimental data.
Acknowledgements
The authors wish to acknowledge
the finnancial support provided by FAPESP, which granted Dr. de Matos a scholarship,
and FINEP and CNPq, for the financial aid for instrumenting the multiphase flow
loop. They all are sponsoring agencies in Brazil.
References
Antal S.P., Lahey Jr, R.T. e Flaherty,
J.E., 1991, Analysis of phase distribution in fully developed laminar bubbly
two-phase flow, Int. J. Multiphase Flow, vol. 17, n. 5, pp. 635-652.&&&&&&&&[  ]Chen, K., 1995, Comparison of different
k-e models for indoor air flow computations, Num. Heat Transfer part B, vol.
28, pp. 353-369&&&&&&&&[  ]Dias S.G., Fran&a F.A. e
Rosa E.S., 2000, Statistical method to calculate local interfacial variables
in two-phase bubbly flows using intrusive crossing probes, Int. Journal of Multiphase
Flow, 26, pp. .&&&&&&&&[  ]Drew, D. and Lahey Jr, R.T., 1979,
Application of general constitutive principles to the derivation of multidimensional
two-phase flow equations, Int. Journal Multiphase Flow, vol. 5, pp. 243-264.&&&&&&&&[  ]Drew, D.A. and Lahey Jr, R.T., 1987,
The virtual mass and lift force on a sphere in rotating and straining inviscid
flow, Int. Journal Multiphase Flow, v.13, n. 1, pp. 113-121.&&&&&&&&[  ]Hibiqui T., Mishima K., Yoneda K.
Fujine S. Tsuruno A. Matsubayashi M., 1994, Visualization of fluid phenomena
using a high frame-rate neutron radiography with a steady thermal neutron beam,
Nuclear Instruments and Methods in Physics Research a vol. 351 pp. 423-436&&&&&&&&[  ]Ishii, M., 1975, Thermo-Fluid Dynamic
Theory of Two-Phase Flow, Eyrolles, Paris - France&&&&&&&&[  ]Keska, J.K. e Fernando, R.D., 1994,
Average physical parameters in an air-water two-phase flow in a small, square-sectioned
channel, J. Fluids Engineering, Vol. 116, pp. 247.&&&&&&&&[  ]Kuo, J.T. and Wallis, G.B., 1988,
Flow of bubbles through nozzles, Int. J. Multiphase Flow, Vol. 14, n. 5, pp.
547.&&&&&&&&[  ]Lahey R.T., 1990, The analysis of
phase separation and phase distribution phenomena using two-fluid models, Nuclear
Engng. & Design, Vol. 122, pp. 17.&&&&&&&&[  ]Lahey, Jr., R.T., Lopez de Bertodano,
M. and Jones, O.C., 1993, Phase distribution in complex geometry conduits, Nuclear
Engineering and Design, vol. 141, pp. 177-201.&&&&&&&&[  ]Lance, M. and Bataille, J., 1991,
Turbulence in the liquid phase of a uniform bubbly air-water flow, J. Fluid
Mechanics, vol. 222, pp. 95-118.&&&&&&&&[  ]Lee, S. J., Turbulence modelling
bubbly two-phase flows, PhD Thesis, Rensselaer Polytechnic Institute, USA, 1988.&&&&&&&&[  ]Lopez de Bertodano, M.L., Lahey,
R.T. and Jones, O.C., 1994, Phase distribuiton in bubbly two-phase flow in vertical
ducts, Int. J. Multiphase Flow, vol. 20, no. 5, pp. 805-818.&&&&&&&&[  ]Lopez de Bertodano, M.L., Lahey,
R.T. and Jones, O.C., 1994, Turbulent bubbly two-phase flow data in a triangular
duct, Nuclear Engineering and Design, vol. 146, pp. 43-52.&&&&&&&&[  ]Lopez de Bertodano, M.L., Lahey,
R.T. and Jones, O.C., 1994c, Development of a k-e
model for bubbly two-phase flow, Trans. ASME, vol. 116, march, pp. 128-134.&&&&&&&&[  ]Mari&, J.L., 1987, Modelling
of the skin friction and heat transfer in turbulent two-component bubbly flows
in pipes, Int. J. Multiphase Flow vol. 113, n? 3, pp 309-325&&&&&&&&[  ]Matos, A., Rosa, E. S., Franca,
F. A., Morandim, M. and Morales, R., 1999, Vertical upward bubbly air-water
flow in a vertical square section channel (in portuguese), Proceedings
of the XV COBEM (CD-ROM), &Aguas de Lind&ia, SP, Brazil.&&&&&&&&[  ]Mishima, K. Hibiki, T. e Nishihara,
H., 1993, Some characteristics of gas-liquid flow in narrow rectangular ducts,
Int. J. Multiphase Flow, Vol. 19, # 1, pp. 115-124.
&&&&&&&&[  ]Mishima, K., Hibiki, T.,
Yoneda, K., Fujine, S., Kanda, K., Nishihara, H. e Chung, N., 1994, Annu. Rep.
Res. Reactor Inst. Kyoto University, Vol. 27, pp. 245-253.&&&&&&&&[  ]Moujaes, S e Dougal R. S., 1987,
Experimental investigation of cocurrent two-phase flow in a vertical rectangular
channel, Canadian J. Chem. Eng. Vol 65, October.&&&&&&&&[  ]Sadatomi M., Sato Y., e Saruwatari
S., 1982, Two-phase flow in vertical noncircular channels, Int. J. Multiphase
Flow, vol. 8, n. 6, pp. 641-655.&&&&&&&&[  ]Shiralkar, B. and Lahey, Jr., R.
T.1972, Diabatic local void fraction measurements in Freon 114 with hot-wire
anemometer, ANS Transactions, 15(2).&&&&&&&&[  ]Spalding, B.D., "The Phoenics
Encyclopaedia", 1994.&&&&&&&&[  ]Stuhmiller, J.H., 1977, The influence
of interfacial pressure forces on the character of two-phase flow model equations,
Int. J. Multiphase Flow, Vol.3, pp 551.&&&&&&&&[  ]Troniewski, L. e Ulbrich, R., 1984,
Two-phase gas-liquid flow in rectangular channels, Chem. Engineering Science,
Vol. 39, # 4, pp. 751-765.&&&&&&&&[  ]Wang S.K., Lee S.J., Jones Jr O.
C., e Lahey Jr R. T., 1987, 3-D turbulence structure and phase distribution
measurements in bubbly two-phase flow, Int. J. Multiphase Flow, vol. 13, n 3,
pp. 327-343.&&&&&&&&[  ]&
Paper accepted September, 2004. Technical Editor:
Aristeu Silveira Neto.}