Group Title: 7th International Conference on Multiphase Flow - ICMF 2010 Proceedings
Title: 2.4.5 - Atomization of Stratified Turbulent Flow
ALL VOLUMES CITATION THUMBNAILS PAGE IMAGE ZOOMABLE
Full Citation
STANDARD VIEW MARC VIEW
Permanent Link: http://ufdc.ufl.edu/UF00102023/00056
 Material Information
Title: 2.4.5 - Atomization of Stratified Turbulent Flow Computational Techniques for Multiphase Flows
Series Title: 7th International Conference on Multiphase Flow - ICMF 2010 Proceedings
Physical Description: Conference Papers
Creator: Smith, E.W.
Roberts, R.M.
McCready, M.J.
Publisher: International Conference on Multiphase Flow (ICMF)
Publication Date: June 4, 2010
 Subjects
Subject: stability
atomization
turbulence
stratified
wavy
Orr-Sommerfeld
k-epsilon
oil-gas pipeline
 Notes
Abstract: A fully computational approach to interfacial stability in stratified turbulent flows is presented. The approach is valid for a wide range of conditions from laminar to highly turbulent (Re = 107). The basestate is modelled using a k-epsilon turbulence model that has been tailored to stratified flow including wave-roughening effects. The stability analysis is performed by solving a modified Orr-Sommerfeld equation with a Chebyshev-Tau method. The determination of atomization and other phase transitions comes from calculating the growth rates of long and short waves as they provide the mechanisms for splashing, droplet shearing, and growth of non-linear waves and slugs. Calculations are in good agreement with data from 1” channels at atmospheric conditions and 2” pipes at 2 atm with the working fluids being water and air. But more than that, the calculations determine the growth rate of detectable waves in experiments with finite residence times.
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: VID00056
Source Institution: University of Florida
Holding Location: University of Florida
Rights Management: All rights reserved by the source institution and holding location.
Resource Identifier: 245-Smith-ICMF2010.pdf

Full Text

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


Atomization of Stratified Turbulent Flow


E.W. Smith*, R.M. Robertst and M. J. McCready*

Department of Chemical and Biomolecular Engineering, University of Notre Dame, Notre Dame, IN 46556, USA
t Chevron Energy Technology Company, Houston, TX 77002, USA
mjm@nd.edu and roberrml@chevron.com

Keywords: stability, atomization, turbulence, stratified, wavy, Orr-Sommerfeld, k-epsilon, oil-gas pipeline




Abstract

A fully computational approach to interfacial stability in stratified turbulent flows is presented. The approach is valid for a
wide range of conditions from laminar to highly turbulent (Re = 107). The basestate is modelled using a k-epsilon turbulence
model that has been tailored to stratified flow including wave-roughening effects. The stability analysis is performed by
solving a modified Orr-Sommerfeld equation with a Chebyshev-Tau method. The determination of atomization and other
phase transitions comes from calculating the growth rates of long and short waves as they provide the mechanisms for
splashing, droplet shearing, and growth of non-linear waves and slugs. Calculations are in good agreement with data from 1"
channels at atmospheric conditions and 2" pipes at 2 atm with the working fluids being water and air. But more than that, the
calculations determine the growth rate of detectable waves in experiments with finite residence times.


Introduction

Multiphase flows occur commonly in the petroleum
production and transportation infrastructure. Despite many
experiment studies and complete knowledge of the
fundamental governing equations, it is still not possible to
predict with confidence pressure drops, slugging and
atomization conditions and solids carrying capacity for
many important operating ranges. This paper discusses
our new computational efforts to predict the onset of
atomization for gas-liquid flows. In addition to altering
the pressure drop flow rate relation, atomized drops can
wet the top of a pipe and hence reduce top-line corrosion.
Our approach is to calculate the basestate (turbulent gas,
laminar or turbulent liquid) for both a smooth surface and
wave-roughened conditions. Stability analysis is done for
both basestates so that secondary transitions, such as the
onset of roll waves, can be determined. Atomization has
been associated with roll wave onset in some previous
studies and atomization can also occur on thin films if
solitary waves are present. Our results consider horizontal,
inclined and declined flows. This work goes beyond
previous studies in some important respects. First it is
recognized that when comparing calculated and measured
basestates for gas-liquid flows, the interfacial roughness
must be considered carefully. Second, the onset of long
waves, since these generally occur at conditions when short
waves are unstable, need to be calculated from this already
rough basestate.

The first question to be addressed is the accurate prediction
of smooth and wavy base-states. A finite difference
method is used to solve a k-epsilon turbulence model for the
gas and liquid flows. This approach will be shown to match
available data for smooth and rough interfaces very well.


The physical phenomenon that occurs in the presence of
short waves is a vertical shift of the maximum velocity
away from the interface. This occurs because the wavy
interface dampens momentum, which can be modeled as an
increase in the eddy viscosity at the interface. These
considerations have not been made in the literature
pertaining to stability and yet they are the key initial step in
the investigation of long wave formation.

A stability analysis was performed using an
Orr-Sommerfeld equation that was modified to include the
basestate velocity and the eddy viscosity. It was found that
using the Chebyshev-Tau method on the modified
Orr-Sommerfeld equation at higher Reynolds numbers and
lower wavenumbers became increasingly difficult to resolve
numerically. To check the validity of the solutions being
obtained, a long-wave analysis was implemented. The
long-wave analysis gives an efficient means of calculating
the neutral stability lines while the Chebyshev-Tau method
provides information about the growth rates of short waves,
which can be correlated to the magnitude of interfacial
roughness giving a completely theoretical prediction of the
onset of atomization. The stability boundary for long
wave predictions has been compared to data and
consistently under-predicts the roll wave transition as shown
in the middle figure. However, this should be the case
because the data are from channels of lengths no more than
10m whereas the stability analysis produces the neutral
stability curve, which corresponds to waves that would
require an infinite length to develop. A contour plot of
growth rates reveals that the roll wave transition matches
roughly with a predicted rate of about 0.1/s which is
nominally the rate necessary for roll waves to form in the
-10m channels where the measurements were taken. Also
in this case, atomization is seen to occur with the presence





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


of roll waves.

Nomenclature


cx Turbulence model constants
c Wave speed (m.s'1)
d Channel height or pipe diameter
f Smoothing turbulence model function
F Froude number
g Gravitational constant (m.s 2)
h Perturbation height (m) (in stability)
h Holdup (in basestate)
m Viscosity ratio
p Perturbation pressure (Pa)
P Pressure (Pa)
r Density ratio
R Reynolds number at interface (Uo.dL.n1)
Re Reynolds number using avg phase velocity
U Basestate velocity (m.s1)
v* Friction velocity (m. s-')
T Weber number
x Position in the flow direction

Greek letters
a Wave number (m'1)
y Perturbation variable
F Total viscosity (Pa.s)
E Turbulent dissipation (m2.s-3)
k Turbulent kinetic energy (m2/s-2)
1T Eddy viscosity (Pa.s)
p Dynamic viscosity (Pa.s)
p Density (kg.m3)
p Stream function (m2.s-')
T Shear stress (N.m 2)


numerical considerations.


Because we are interested in the stability of the turbulent,
potentially roughened basestate, we need an approach for
which stability can be implemented. Methods such as
Direct Numerical Simulation or Large Eddy Simulation are
not only extremely time consuming, but it is doubtful that
information from these calculations could readily transfer to
a stability analysis. Alternatively, highly empirical
methods such as zero equation models or mixing length
models are not flexible enough or lack the empirical values
to use with high pressure oil-gas flows. Two equation
models were chosen for their application to a wide range of
conditions and also in the flexibility afforded by them to
adjust for the effects of roughness and turbulent effects at
low Reynolds numbers.

The k-epsilon, two-equation models are a way of
determining the eddy viscosity, which is part of the
Boussinesq approximation. By modelling the Reynolds
stress as a second viscosity, the momentum equation
simplifies to


dy(a )d ) = Pdx
ay ~ d d jP Slr)


where p, is p, = pjc,, k2/E .


The modeling


parameters c, and f, depend on version of the turbulence
model and Z is the isotropic part of the turbulent
dissipation. The variables, k and E are determined by
the continuum equations


d Ty dk y Pj (E+ D)
kdy k dy


Subscripts
i Interface
j Phase index (L or G)
k Turbulent kinetic energy
L Liquid phase
G Gas phase
e Turbulent dissipation
p Eddy viscosity equation

Superscript
+ Dimensionless friction variable


Basestate


Theory
Here we will describe the approach used to model wavy
flows and give an overview of the numerical methods. A
stability analysis is only a tool that provides information
regarding a particular basestate. Thus without an accurate
model of the base flow there is simply no relevance of
hydrodynamic stability results to a desired physical situation.
To model a stratified turbulent flow, a single phase
turbulence model must be adjusted for the roughening effect
of interfacial waves and also increased turbulence in the
liquid at Reynolds numbers well below the empirical
transition (Re=2300). In addition, there are several


d + LT\ dE\ Lt2 dU
=. r C ) cE2 f2 2 C f r PT L )
dy \ E dy )k k (dy)


PjE,


where again cx, fx, E, and D vary depending on the
model. The various models of Chien (1982), Dutoya and
Michard (1981), Hoffman (1975), Lam and Bremhorst
(1981), Launder and Sharma (1974), and Reynolds (1976)
were tested. The version chosen for this work was based on
the Launder and Sharma model, but with f, taken from
Chien.

Results
The maximum velocity of a fluid will move away from a
wavy interface or a rough wall because there is more
resistance to motion in this region. This is modeled simply
as a non-homogeneous boundary condition for the turbulent
kinetic energy. The value of the condition is determined
from matching data, which includes profile data when
available or by matching flow rates and pressure drops. It
was found that in general k+ = k/v*2 = 4 gave a good fit,
where v* is the friction velocity. An example of the
success of this approach is shown in figure 1.


/dU \2
PT -
(dy)














0.6

0.4




0
0.2


0 2 4 6 8
Velocity [mrs]

Figure 1: Velocity profile from Cohen and Hanratty (1968)
versus the model of wavy flow

Further evidence of the viability of the two-equation model
is given in table 1 where the pressure drop data calculations
are compared to data from Bruno (1988). While the
condition, k = 4 is useful, it is not expected to be a
completely general result. It could probably be improved
by using the growth rate of the smooth basestate to estimate
a wave height.

ReG ReL Pressure Drop Pressure Drop error
Data [Pa/m] Calculation

9000 720 27.4 28.2 3%
7670 720 20.5 19.3 -6%
6100 720 14.7 14.9 1%
11000 745 26.2 22.8 -13%
10050 745 23.6 20.5 -13%
9000 745 19.0 17.7 -7%
7670 745 13.9 12.9 -7%
6100 745 9.7 9.4 -3%
4750 745 6.4 6.4 0%
11000 255 18.7 15.6 -17%
10050 255 16.3 14.4 -12%
9000 255 12.6 11.5 -9%
7670 255 9.1 8.8 -3%
6100 255 6.5 6.2 -4%
4750 255 4.3 4.1 -4%
Table 1: Basestate comparison with data from 1" air-water
channel data

Before discussing stability, it is worth reiterating the
importance of the basestate. One goal of this work is to
correctly predict flow pattern transitions. This involves
correctly calculating the stabilizing and destabilizing
mechanisms that include gravity, shear stress, pressure, and
surface tension. If the basestate velocity or eddy viscosity
profiles are not correct then the relative magnitudes of these
mechanisms will give erroneous transition lines. The
transitions of interest are that of stratified to wavy, the onset
of entrainment, the initiation of slug flow, and maximum
entrainment. A turbulent, linear, stability analysis provides


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

crucial information for all of these transitions. Obviously,
non-linear techniques are required for a study of slug flow,
but again such an undertaking is not possible without a
correct linear analysis. Figure 2 illustrates this by showing
how much the wave speeds of long waves depend on the
details of the linear analysis and the underlying basestate.
This has not been available and is why this work poses great
opportunities in producing far more accurate and
fundamentally sound non-linear analyses. Ultimately, it is
the necessary first step in reassessing flow pattern
boundaries, a subject that has seen little improvement since
Taitel and Dukler over 30 years ago.


0.4


0.3


I0.2
-Turbulent
0.1 -Laminar
-Inviscid
0.0
0,01 0.1 1 10 100 1000
Wavenumber Il/m]
Figure 2: Comparing wave speeds of simple laminar and
inviscid flows versus a complete turbulence model and
stability analysis

Linear Stability Analysis

Theory
The Orr-Sommerfeld equation is an expansion of the
laminar terms of the Navier-Stokes equations. Because the
Boussinesq approximation is employed, the expansion can
be repeated with the eddy viscosity term. First, the
velocity is recast by the Reynolds decomposition,
u(x, y, t) = U(y) + yu(x, y, t) and v(x,y,t) = V(y) +
yv(x, y, t), where y is the expansion parameter and is
much smaller than one. The perturbation velocity is
assumed to have the form of a traveling wave, un(x, y, t) =
u (y) exp(ia (x c t)). Use of a stream function means
the perturbation velocities can be written as


M(x,y, t) = 0'(y) exp(ia(x ct))
v(x, t) = ia 0(y) exp(ia(x ct)).


The Boussinesq approximation remains in the momentum
equations, but the eddy viscosity is not decomposed into
basestate and perturbation viscosities. The reasons for this
are twofold, the first is that this term is presumably small
even compared to the other perturbation values (the correct
phase shift of the shear stress confirms this). Also, there is
little reason to believe that the perturbation of a modeling
convenience such as the eddy viscosity accurate as it may
be for the basestate would hold any physical relevance in
the stability analysis. The pressure is also represented as a
sum of a basestate and normal, traveling-wave, perturbation,
p = P + yp(y) exp(ia(x c t)) Upon substituting
these expansions the order y' terms for the horizontal and
vertical momentum equations are










dy
dy (r(y)i "(y)) 2r' (y)
rk
= aiR R (p(y) c Y'(y) P(y)U(y)
mk
+ )'(y)U(y))

and

a(a22(y)P(y) '(y)'(y) f(y)v"(y))
= i r (a2C 2 (y) p(y)
mk
a2 U(y)P(y)),

respectively. Differentiating the horizontal momentum
equation with respect to y, multiplying the vertical
momentum equation by a, and adding the two equations,
the pressure perturbation can be eliminated giving

d2 d
d ((y)"(y)) 22 d~ (r(y)('(y)) + a4r(y)(y)
dy2 dy
-iRk ((U(y)- c)(a2 (y)
mk
"(y)) + U"(y)(y)),

which is a modified Orr-Sommerfeld equation. Unlike the
standard form, it includes a basestate eddy viscosity as well
as the basestate velocity profile.

In addition to an expansion of the differential equations, the
boundary conditions must also be reevaluated to include the
eddy viscosity. The no slip and zero flux conditions at the
walls remain homogeneous:

4=0
(' =0.

Disturbances in the flow produce perturbations at the
interface. Therefore, the boundary conditions do not occur
precisely at y = 1, but on an undulating interface of the
form h(x, t) = yh exp(ia(x c t), where h is a complex
scalar. The displaced velocities are written as a Taylor
series expansion with respect to h around y = 1 giving

du(1)
u(1 + h) = u(1) + h d-
dy
dv(l)
v(1 + h) = +h
dy

Substituting the normal wave expansion into the velocity
matching condition gives

PL = (G
P'L(y) + UL(y)h = b (y) + U'(y)h.

Matching the shear stress requires both horizontal and
vertical velocity components of the Cauchy stress tensor
r = F (du/dy + dv/dx) Expanding and equating this
term across the interface gives

FrL(L' h U (' + aPL) = m2F (PG + h UL' + a2PL).


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


Finally, the normal stress is required to match at the
interface. In addition to the velocity gradients, a
perturbation pressure, the capillary pressure produced from
surface tension, and hydrostatic pressure created by an
uneven interface must be included. The normal stress
condition is


(FL dvL
2 L-y
dy


m y R (PL PG)
dy


dh2 (, t)
S-TRd2
dx2


R F (1 r2)h,


where T is the Weber number, and F is the Froude
number. The pressure for phase j is

P= i- (a2p + 1 (3)
ryR
+ a(U + c ().

Substituting the normal wave expansion, the pressure, and
collecting y' terms gives the normal stress condition

3ia2
-T a3h aF h (1 r2) + ( ; m2FG')
+ a(#PU' Uj(P rT2ct 2 r222
1 d d
Rd (L L')- m2 y (FG )
= -ca(P' r2 9).

The system is still underdetermined. There are 8 boundary
conditions, 2 fourth order equations, and an unspecified
interfacial perturbation variable h The kinematic
condition provides the necessary condition for solvability

h = c h.

In practice, one of the boundary conditions is dropped by
solving for h using one of the above conditions. While
much of the literature uses the kinematic condition to
determine h, this leads to an eigenvalue in the denominator
and thus the matrix formulation does not work. Therefore,
the value of h determined by the continuity of tangential
velocity is substituted into all the boundary conditions. It
is possible to keep h as a variable in the vector of
unknowns, and use all 9 boundary conditions. In this case
the vector would include the value of the series coefficients
followed by the perturbation height.

Results
The solution of the eigenvalue problem is to factor the
modified Orr-Sommerfeld equation into two second order
differential equations and apply the Chebyshev-Tau Spectral
method. The application of this method is not
straightforward. The difficulty is using enough ( Ii.bI sheIc\
polynomials to resolve the basestate profiles and their
derivatives, but few enough that numerical error remains
insignificant. For moderate flows this condition is easily
met. Figure 3 shows clearly defined growth rate and wave
speed curves for the interfacial modes of a moderate flow of
water and air in a 5" channel.
























10 102
wavenumber [1/rnm]


100 102
wavenumber [1/m]
Figure 3: Growth rate and wave speeds of interfacial modes
at ReG=31,000 and ReL=3900. The interfacial modes
correspond to the downstream (solid) and upstream (dashed)
traveling waves

Increasing the inertia necessitates a finer grid with more
points in the basestate, which in turn requires more modes in
the stability problem. The consequence is a narrow or
non-existent range of mode-node combinations that produce
realistic profiles. To broaden the space between poor mode
expansions and machine precision, several numerical tricks
were implemented as described in the thesis of Smith (2010).
The success of these numerical techniques is shown in
figure 4.

When data are unavailable, there are other implicit means of
checking the validity of results. For instance, it is not
sensible to say that upstream traveling long waves will grow
while downstream traveling waves are stable, which is what
the upper plot of figure 4 indicates. For some conditions
well above stability a linear analysis may reveal that
upstream waves are unstable, but this is only valid if the
downstream modes are unstable first. Intuition can serve
as a check in other ways. An example is that the speed
curve should reach a minimum where the growth rate is
highest, this is a well-known energy trade-off between
growing versus propagating the wave. Additionally, the
interfacial modes should be mirror images of one another, a
mathematical and physical truism of wave behavior. The
speeds should also increase where the restoring force of
surface tension (large wavenumbers) and gravity (small


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

wavenumbers) are applicable. The bottom curves of figure
4 pass all of these tests and are easily defended as valid
growth rates and wave speeds.


U)
E0.5



,-0.5

-1



2


E0.5

_0

-.5
-1


100 102
wavenumber [1/m]


100 10/
wavenumber [1/rm]


Figure 4: For highly turbulent flows (ReG=3.2x105,
ReL=4000) the solution becomes exceedingly difficult. By
judiciously implementing a number of numerical techniques
the mode interactions in (a) can be eliminated (b).

Having developed a stability method that is at the very
minimum defensible (the accuracy will be defended in the
next section), this method can be used to probe parameter
space and determine where instabilities might arise due not
only to flow rates but flow parameters such as viscosity
ratio, density ratio, holdup, pipe diameter, incline, etc. An
example of this type of analysis is given in figure 5. In the
first plot the increase of pressure, while maintaining the


'N"-

/









same Gas Reynolds numbers, decreases stability and
ultimately leads to the formation of long waves. Similarly,
changing the flow from a decline to an incline produces a
long wave transition.


P=1atm
- P=2atm
P=1.5atm


0-1

-0.212
102


S0.3 r 9=+0.03
0.2 -

0.1

D o0 /

-0.1

-0.2
10 10 102
a
Figure 5: Long wave transitions produced by increasing a) the
system pressure and b) the inclination


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

was calculated and the resulting growth rate contours are
compared with Bruno's long wave transition in figure 6.


10000


5000


500 1000
ReL


1500 2000


Figure 6: Bruno's experimental long wave transition
(dashed line) as determined by a spectral analysis of wave
tracing is compared to the growth rate contours of long
waves as determined by a stability analysis

We can use growth rate contours to find neutral stability for
highly inertial flows. A program has been developed that
handles such conditions and produces contour plots. Some
data are available that determine the onset of atomization by
visually confirming the presence of droplets in air. Since
these experiments were done in relatively short pipes it
should be expected that growth rates corresponding to onset
would be higher than that of Bruno's spectral analysis. A
plot of growth rate contours and the onset of atomization for
a water-air system in a 2" by 20' pipe is given in figure 7.
The kinks in the lines result from the tradeoff between
precision and speed between automation and a laborious
adjustment of numerical parameters. If need be, trend
lines can be fitted to smooth contours which is shown in
figure 8.


Phase Transition and Atomization


0.1


A motivation of this work is to study the onset of
atomization. A mechanism for this is the shearing of liquid
from the crest of roll waves. Therefore a prediction of
long wave instability should lead to waves that will correlate
with the presence of atomization. Here we compare the
stability calculations of long waves with experimental
values for the transition. The difficulty in measuring the
existence of roll waves is that laboratory experiments are
performed in channels or pipes no longer than about 10 m,
and thus detectable waves are limited by growth rates and
finite residence times. Furthermore, outside of some
heuristics, there is no quantitative method of determining
the growth rate of these realizable waves. A stability
analysis provides a quantitative context for experimental
observations.

At conditions or flow distances when roll waves are not
visible, they can still be detected by a spectral analysis of
the wave tracings. Bruno and McCready (1988) performed
this analysis for water-air experiments in a 1" inch channel.
These data are ideal for testing the stability analysis because
the flow rates are moderate and the analysis provides the
best experimental estimate of neutral stability. The
stability of long waves at various points in ReL-ReG space


0.01

0.001

0.0001


Use

Figure 7: Growth rate contours of a highly inertial flow and
the experimental onset of Mantilla (2008)


~6)

~ bC


.5 3

.2 \
0.1




.I . . . I













01


001


0.001


0.0001


- \
\ 3.0
0.5 1









).1 1 10 100


USG

Figure 8: Fitted growth rate contours

The stability analysis also includes the growth rates of short
waves. One reason to doubt the applicability of the long
wave mechanism to all systems is the evidence in data that
surface tension has an effect on the transition line. Since
surface tension only appears in one boundary condition and
at order a3 where a is the wavenumber, it is likely that
long wave instability is not the only criterion for
atomization. Therefore, the short wave growth rate
contours can be compared to data with the purpose of
finding a growth rate threshold that corresponds to the
observed onset of atomization. An example of this
analysis is given in figure 9. The short wave analysis has
only recently been considered, but it is presented here as an
example of information that is only available with a
comprehensive approach to turbulent stability analysis.


I


0.1


0.01


0.001


0.0001


Figure 9: Short wave growth rates
experimental onset line (o)


compared to the


Another possible application of this code is determining the
maximum entrainment. The onset of entrainment is the
point where the first liquid drop becomes airborne. At
maximum entrainment there is no net increase in the liquid
content of the gas phase. At this point the interface is once
again stable with respect to long waves. The progression
from onset to maximum entrainment can be posed as a
stability problem applied to a mass balance. In other
words, given two flow rates the percent entrained
determines the gas density, the liquid height, and the
stability. Simply moving from zero to 100% entrained will
produce a line such as the one found in figure 10, which
determines the location of the onset and maximum
entrainment for those flow rates.


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


x 10





.o0
W- 2 %00
0
ct' 0


1 0 Unstable
0 Stable Below Onset
O Stable Max Entrainment

0 50 100 150 200 250 300 350
Re,
Figure 10: Entrainment line for a liquid mass flow rate of
0.53kg/min and a gas mass flow rate of 0.10kg/min. The
50 points span entrainment values from 0 to 99%

Summary

This paper briefly describes a new stability procedure for
turbulent gas-liquid flows. A k-epsilon model is shown to
provide good predictions of wavy stratified flows. The
onset of atomization is predicted by the instability of long
waves that occur on roughened two-layer basestate. The
criterion is not neutral stability but a finite growth rate
sufficient to allow such waves to grow in finite length
systems. The conditions of maximum entrainment are
predicted when the longwaves again become stable because
the liquid layer is getting thinner as the entrainment gets
larger. This boundary is calculated by reducing the liquid
flow rate and adjusting the gas flow density and viscosity to
account for the entrained drops.

Acknowledgements

The authors would like to thank the Chevron Energy
Technology Company for their funding.

References


Bruno, A Study of Interfacial Waves in Gas-Liquid Flows,
o University of Notre Dame (1988)

Bruno K., McCready M.J., "Origin of roll waves in
horizontal gas-liquid flows" AIChE J. 34, 9, 1431-1440
(1988)

Chien K.Y, Predictions of Channel and Boundary-Layer
Flows with a Low-Reynolds-Number Turbulence Model,
AIAA Journal, 20, 33-38, (1982)

Cohen L.S., Hanratty T.J., "Effect of Waves at a Gas-Liquid
Interface on a Turbulent Air Flow" J. Fluid Mech. 31, 3,
467-479 (1968)

Dutoya D., Michard P., "A Program for Calculating
Boundary Layers along Compressor and Turbine Blades",
Numerical methods in Heat Transfer, John Wiley & Sons,
New York (1981)

Hoffman G.H., "Improved Form of the Low-Reynolds
Number k-epsilon Turbulence Model", Physics ofFluids, 18,
309-312 (1975)


S20 40


0 '
0
5 "





. . . . . -. . .






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

Lam C.K.G., Bremhorst K.A., "Modified Form of the
k-epsilon-Model for Predicting Wall Turbulence", Journal
of Fluids Engineering, 103, 456-460 (1981)
Launder B.E., Sharma B.I., "Application of the
Energy-Dissipation Model of Turbulence to the Calculation
of Flow Near a Spinning Disc", Letters in Heat and Mass
Transfer, 1, 131-138 (1974)
Mantilla I., "Mechanistic modeling of liquid entrainment in
gas in horizontal pipes". University of Tulsa (2008)
Reynolds W.C., "Computation of Turbulent Flows", Annual
Review ofFluid Mechanics, 8, 183-208 (1976)
Smith E.W., Thesis, University of Notre Dame (2010)
Unpublished
Taitel Y., Dukler A.E., "A Model for Predicting Flow
Regime Transitions in Horizontal and Near Horizontal
Gas-Liquid Flow", AIChE J, 22, 1, 47-55 (1976)




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 - - mvs