Revista Mexicana

de Ingeniera Qumica
Academia Mexicana de Investigaci
on y Docencia en Ingeniera Qumica, A.C.

Revista Mexicana de Ingeniera Qumica

Volumen 11, N
umero 2, Agosto 2012

Vol. 11,CONTENIDO
No. 2 (2012) 269-278

ISSN 1665-2738

Volumen 8, nmero 3, 2009 / Volume 8, number 3, 2009

RISING VELOCITY FOR SINGLE BUBBLES IN PURE LIQUIDS

VELOCIDAD DE ASCENSO DE BURBUJAS AISLADAS EN LIQUIDOS

PUROS
213 Derivation and application
of the Stefan-Maxwell
equations 2 and A. Soria1
S. Baz-Rodr
guez1 , A. Aguilar-Corona
1

Departamento
de Ingenier
a dede
Procesos
e Hidrde
aulica,
Universidad Autonoma Metropolitana - Iztapalapa, San
(Desarrollo
y aplicacin
las ecuaciones
Stefan-Maxwell)
Rafael
Atlixco
186,
Col.
Vicentina,
Mexico, D. F. CP. 09340, Mexico.
Stephen Whitaker
2
Facultad de Ingeniera Mecanica, Universidad Michoacana de San Nicolas de Hidalgo, Francisco J. Mujica s/n,
Ciudad Universitaria, CP. 58040, Morelia, Michoacan, Mexico.
Biotecnologa / Biotechnology
Received October 9, 2011; Accepted March 10, 2012.
245 Modelado de la biodegradacin en biorreactores de lodos de hidrocarburos totales del petrleo

Abstract
An equation
to predict theenterminal
velocity of single bubbles in stagnant Newtonian liquids is proposed. The
intemperizados
suelos y rise
sedimentos
formulation combines a force balance obtained from the boundary layer theory for non-distorted bubbles and an
analytic equation
coming modeling
from a mechanic
balance.
including
geometric
parameters,
(Biodegradation
of sludge energy
bioreactors
of totalWithout
petroleum
hydrocarbons
weathering
in soil which are
difficult to and
assess,
it
is
assumed
that
the
weighting
of
dominant
forces
is
enough
to
adequately
predict
the terminal
sediments)
velocity in both the intermediate and inertial motion regimes. The proposed equation shows good agreement with
S.A.data
Medina-Moreno,
Huerta-Ochoa,
C.A. Lucho-Constantino,
L. Aguilera-Vzquez,
Jimnezexperimental
from bubblesS.rising
in pure liquids.
Moreover, for bubbles
rising in cleanA. water,
the effect of
helical trajectories
from experimental data trends and included in the formulation as a correction
Gonzlez was
y M. estimated
Gutirrez-Rojas
factor for
terminal velocity.
259the
Crecimiento,
sobrevivencia y adaptacin de Bifidobacterium infantis a condiciones cidas
Keywords: single bubble motion, terminal velocity equation, oscillatory bubble path, pure liquids.
(Growth, survival and adaptation of Bifidobacterium infantis to acidic conditions)

Resumen
Mayorga-Reyes,
P. Bustamante-Camilo,
A. Gutirrez-Nava,
E. Barranco-Florido
A. AzaolaSe propusoL. una
ecuacion para
predecir la velocidad
terminal del ascenso
de burbujasy aisladas
en lquidos
Newtonianos.
La formulacion combina un balance de fuerzas obtenido de la teora de capa lmite para burbujas
Espinosa
esfericas
una ecuaci
on anal
proveniente
un balance
de energ
mecanica. Sin cerevisiae
la inclusioin
n de
265con
Statistical
approach
totica
optimization
of de
ethanol
fermentation
byaSaccharomyces
theparametros
geometricos que son difciles de determinar, se asume que la ponderacion de las fuerzas dominantes es suficiente
para predecir
adecuadamente
la velocidad
presence
of Valfor zeolite
NaA terminal en los regmenes de movimiento intermedio e inercial. Las
predicciones de la ecuacion propuesta presentan buen ajuste con datos experimentales para burbujas ascendiendo en
(Optimizacin
de laascendiendo
fermentacinenetanlica
de Saccharomyces
cerevisiae
presencia de
lquidos puros.
Ademas, estadstica
para burbujas
agua pura,
se estimo el efecto
de lasen
trayectorias
helicoidales
a partir de zeolita
la tendencia
de
datos
experimentales
y
se
incluy
o
en
la
formulaci
o
n
como
un
factor
de
correcci
on de la
Valfor zeolite NaA)
velocidad terminal.
G. Inei-Shizukawa, H. A. Velasco-Bedrn, G. F. Gutirrez-Lpez and H. Hernndez-Snchez

Palabras clave: movimiento de burbujas aisladas, ecuacion de velocidad terminal, trayectoria oscilatoria de
burbujas, lquidos puros.
Ingeniera de procesos / Process engineering

271 Localizacin de una planta industrial: Revisin crtica y adecuacin de los criterios empleados en
Introduction
velocity defect caused by the wakes (Marks, 1973). If
esta decisin

one is dealing with bubble columns, the slip velocity

between
can be expressed in terms of a single
(Plant of
siteaselection:
Critical
review
adequation criteria
usedphases
in this decision)
Motion analysis
single bubble
rising
in and
quiescent
bubble
terminal
velocity and the gas volume fraction
J.R. Medina,
R.L. Romero
y G.A. Prez
liquid is useful
to understand
and describe
gas-liquid
bubbly flows. Let us consider a train of bubbles rising
unconfined through still liquid. Its rise velocity can
be approximated by that of a single bubble plus the

(Shah et al., 1982; Ruzicka et al., 2001). These facts

emphasize the importance of having equations that
predict the terminal rise velocity of single bubbles, sin-

Corresponding author. E-mail: [email protected]

Tel.: +52 55 58044648 Ext. 237

Publicado por la Academia Mexicana de Investigacion y Docencia en Ingeniera Qumica A.C.

269

Baz-Rodrguez et al./ Revista Mexicana de Ingeniera Qumica Vol. 11, No. 2 (2012) 269-278
ce it is a fundamental parameter in gas-liquid flows.
The shape and terminal velocity of isolated
bubbles vary as the equivalent diameter increases,
according to three regimes (Tomiyama et al., 2002):
1) dominant viscosity, 2) intermediate region, where
viscous, surface tension and inertial effects should
be taken into account, and 3) dominant inertia. Out
of the three regimes, the intermediate one presents
major difficulties since the force interactions and
their effect on the terminal velocity is not yet well
understood.
Moreover, the terminal velocity is
particularly sensitive to the presence of contaminants
and even to the initial deformation induced by
injection (Tomiyama et al., 2002; Celata et al., 2007).
The intermediate regime spans an equivalent
diameter range from 0.07 to 0.65 cm for air bubbles
ascending through pure water (Mendelson, 1967;
Loth, 2010). This corresponds to particle Reynolds
numbers (Re) between 80 and 1500 for normal
pressure and temperature conditions, approximately.
This working interval is usual in many industrial
systems. When Re > 200, bubble shapes are oblate
ellipsoidal with fore-aft symmetry. This symmetry
gradually brakes as Re increases, until the bubble
motion attains the inertial spherical cap regime.
Besides, for Re 650-685, the onset of path instability
occurs, and a transition from rectilinear to zig-zag
or helical trajectories is observed (Duineveld, 1995;
Sanada et al., 2008; Veldhuis et al., 2008).
There exist three basic approaches to establish the
terminal velocity (Kulkarni and Joshi, 2005):
Force balance. The terminal velocity equations
result from a force balance between drag and
buoyancy. The most representative solutions are
obtained from drag in creeping and potential
flows. This approximation can be adequate to
predict the rise velocity of small bubbles when
viscosity effects are still dominant.
Dimensional analysis. Dimensionless groups
are determined from the leading variables that
govern the bubble motion. A functional relation
is proposed for these groups, and the adjustable
parameters are fitted from experimental or
numerical data (Rodrigue, 2001).
Wave analogy. Interfacial disturbances are
assumed for bubbles whose dynamics is similar
to propagating waves in an ideal fluid, and
then the terminal velocity is estimated as a
function of the equivalent diameter and the fluid
properties (Mendelson, 1967). Even if some
270

doubt has been cast on the original approach

because of its lack of physical basis (Lehrer,
1976; Fan and Tsuchiya, 1990; Bozzano and
Dente, 2001), the wave analogy approach can
span some later investigations, which obtained
equations with the same basic form (Lehrer,
1976; Zudin, 1995). Predictions using this
approach fit the main trend after the bubble
size corresponding to the terminal velocity local
maximum is attained (0.14-0.18 cm for air in
water, see Fig. 1).
Other formulations for bubble terminal velocities have
been proposed from some detailed bubble shape and
external flow considerations (Moore, 1965; Bozzano
and Dente, 2001; Tomiyama et al., 2002). However,
these approaches often require additional data such
as aspect ratio (between minor and major axes)
and (hemispherical) bubble distortion factor. These
parameters are difficult to obtain. Some investigations
report the aspect ratio from experimental results on
pure and contaminated liquids as a function of the
Eotvos number (Eo) (Wellek et al., 1966; Okawa et
al., 2003; Sanada et al., 2008). Sometimes the Weber
number (We) is used instead (Moore, 1965; Wellek et
al., 1966), which seems more appropiate to establish
functional relations with this geometrical parameter
(Celata et al., 2007).
Nevertheless, expressing
the aspect ratio in terms of We involves its own
dependency on the terminal velocity. Besides, for
pure liquids it is difficult to estimate aspect ratios
or distortion factors even experimentally (Tomiyama
et al., 2002), in contrast with the case of liquids
contaminated with surfactants, which act by damping
bubble shape oscillations (Celata et al., 2007).
The inclusion of drag coefficients with non-linear
Re dependency, or aspect ratio that depends on We
in any terminal velocity formulation implies solving
non-linear equations, often with fractional exponents.
This is inconvenient when the terminal velocity is not
the main objective of an engineering calculation. For
this reason, a simple equation to predict the isolated
bubble terminal velocity is proposed in the present
work. If geometrical considerations are avoided, we
believe that the proper weighting of the dominant
forces is enough to assess this important parameter for
the intermediate and inertial regimes.
The terminal velocities reported in literature
correspond to the axial component of the bubble
center of mass velocity vector. However, for the
case of air rising through water, it is important to
take into account that rectilinear trajectories become

www.rmiq.org

Baz-Rodrguez et al./ Revista Mexicana de Ingeniera Qumica Vol. 11, No. 2 (2012) 269-278
unstable above some critical Re. For this reason, and
considering a helical trajectory, a correction factor was
obtained as a function of the mean slope of the helix.
This procedure improves the fit to experimental data
in the equivalent diameter range that correspond to the
oscillatory motion of the bubble.

numbers smaller than unity (Tomiyama et al., 2002)

and Mendelson equation seems not to have an
adequate basis to justify its deduction (Lehrer, 1976;
Bozzano and Dente, 2001).
When viscous effects are losing their influence, the
terminal velocity for a bubble in potential flow can be
used (Levich, 1962)

Terminal velocity formulation

VT pot =

Considering the intermediate region, inertial effects

are combined with viscosity and surface tension
effects. This suggests that an equation for the bubble
terminal velocity, which is valid when viscous effects
are still important, can be combined with another one
that is valid when surface tension effects are significant
(Fan and Tsuchiya, 1990; Jamialahmadi et al., 1994).
In the former case it is reasonable to propose an
equation developed from the force balance approach.
In the later it would be sensible to use a wave analogy
approach. Proceeding this way, an equation can be
obtained, which is valid for a wide range of equivalent
bubble diameters.
The proposed formulation for the terminal velocity
is given here by the following combination of the
corresponding effects.
VT = q

1
1
VT2 1

1
VT2 2

(1)

where VT 1 is the rise velocity when viscous effects are

still important and VT 2 is the corresponding velocity
when surface tension effects are significant. The form
of Eq. (1) was suggested by Jamialahmadi et al.
(1994) and remains in use (Kulkarni and Joshi, 2005;
Cai et al., 2010; Duangsuwan et al., 2011); also, a
similar form but with fitting parameters was previously
proposed by Fan and Tsuchiya (1990). When viscous
effects are important and the bubble diameter is small,
Eq. (1) ensures that the contribution of VT 2 is
neglected. This is important because the terminal
velocity equations from wave analogy yield values
too high below the intermediate velocity maximum,
and they are preferably valid until significant surface
tension effects appear. The opposite occurs with VT 1 at
larger bubble diameter. In the intermediate region, Eq.
(1) is just the terminal rising velocity coupled between
the two velocities VT 1 and VT 2 .
Jamialahmadi et al. (1994) proposed the equation
developed by Hadamard (1911) and that of Mendelson
(1967) for VT 1 and VT 2 , respectively. Nevertheless,
Hadamard equation is only valid for Reynolds

1 gde2
36 L

(2)

where is the density difference between the liquid

and gaseous phases, de is the diameter of a spherical
bubble of equivalent volume, L is the dynamic
viscosity of the liquid and g is the acceleration
of gravity. Contrary to Hadamard solution, which
overestimates the terminal velocity for moderate Re,
the potential solution underestimates them. Using
the governing equations of motion for the boundary
layer on a spherical bubble, Moore (1963) obtained
a correction for the potential solution of drag force
which fits particularly well for moderate regime
motion (50 < Re < 200), before the surface tension
and inertial effects start to cause deviations from the
spherical bubble shape. Consequently, the boundary
layer solution must be more appropriate when the
bubble rise at moderate Re.
Buoyancy and drag as expressed by Moore (1963)
gives the following equilibrium:

2.21
(3)
6L de VT 1 1 1/2 = de3 g
6
Re
1

where

de VT 1 L
(4)
L
Rearranging Eq. (3), the solution for the terminal
velocity for potential flow appears
Re1 =

V 3/2
1/2 T 1

VT1/2
pot (gde )

VT1/2
pot
(gde )1/2

VT1/2
1 0.36833 = 0 (5)

Expanding a Taylor series around the potential flow

solution, and truncating to the second order term, an
equation for the terminal velocity is obtained. The
above allows obtaining an approximate solution for the
bubble velocity avoiding numerical solving of Eq. (5).
The resulting equation is as follows
#1/2
"
(gde )1/2
(6)
VT 1 = VT pot 1 + 0.73667
VT pot
Eq. (6) is valid for spherical isolated bubbles rising
under the effect of a hydrodynamic boundary layer and
buoyancy force.

www.rmiq.org

271

Baz-Rodrguez et al./ Revista Mexicana de Ingeniera Qumica Vol. 11, No. 2 (2012) 269-278
When surface tension effects are important, Lehrer
(1976) states a mechanical energy balance generated
by bubble displacement in the liquid:

Eq. (8) can predict the terminal velocity when

the dominant effects come from surface tension and
inertia. According to the combination proposed by Eq.
(1) it is possible to calculate the terminal velocity for
a large range of bubble sizes in pure liquids for the
intermediate and inertial regimes.

Results and discussion

Predictions from Eqs. (6) and (8), combined according

to Eq. (1), were compared with experimental data
from Haberman and Morton (1953), Peebles and
Garber (1953), Okazaki (1964), Aybers and Tapucu
(1969), Duineveld (1995), Blandn-Arrieta (1997),
Leifer et al. (2000), Okawa et al. (2003), Talaia (2007)
and Sanada et al. (2008) for gas bubbles ascending
through various pure liquids, including distilled and
ultra-purified water.
Fig. 1 shows the comparison of the experimental
data for distilled water case in a wide equivalent
diameter range with the terminal velocity predicted
from equations. The proposed formulation in this
work properly predicts the typical behavior of the
terminal velocity along the whole range. Predictions
of three other formulations are also presented. Two
of them propose a weighting similar to that of Eq.
(1) (Fan and Tsuchiya, 1990; Jamialahmadi et al.,
1994), and the other one was obtained by dimensional
analysis (Rodrigue, 2001).
272

6 0

5 0

3 0

(c m /s )

4 0

1 31
1
d L 2 = de2 + de3 gde
(7)
6 e2
6
where is the surface tension and is the liquid
velocity when it is displaced by the bubble passing by.
The right hand side of Eq. (7) is the potential
energy increase of the displaced liquid due to surface
tension and buoyancy during the bubble motion over a
distance de . This potential energy is then transformed
in kinetic energy [left hand side of Eq. (7)] as the
bubble moves, and finally it is dissipated by the wake.
The liquid motion is considered uniformly accelerated.
Assuming that the necessary time for the potential
energy to be transformed into kinetic energy is that
corresponding to a displacement de of the bubble
moving at a steady velocity, the following equation is
obtained (Lehrer, 1976).
!1/2
gde
3
+
(8)
VT 2 =
L de
2L

Haberman and Morton (1953)

Haberman and Morton (1953)
Peebles and Garber (1953)
Okazaki (1964)
Aybers and Tapucu (1969)
Blandn-Arrieta (1997)
Leifer et al. (2000)
Okawa et al. (2003)
Talaia (2007)
Eq. (1) with Eqs. (6) and (8)
Fan and Tsuchiya (1990)
Jamialahmadi et al. (1994)
Rodrigue (2001)

2 0

1 0

0
0 .1

d
e

(c m )

Fig. 1. Comparison between predictions of equations

and experimental data for bubble terminal velocity
through distilled water (7.4 < Re < 32500; 0.002 <
We < 245; 0.01 < Eo < 475) at 20-28 C.
Markers correspond to experimental data, and lines to
predictions from the equations.
The prediction of terminal velocity by
Jamialahmadi et al. (1994) equation fails before the
local maximum, and does not fit the experimental data
at this maximum. This happens because the solution
proposed by Hadamard -and used by Jamialahmadi et
al. (1994) for dominant viscous effects-, overestimates
the terminal velocity value for Re larger than one.
A similar situation occurs with Fan and Tsuchiya
(1990) equation, which is based on a parameter
fit. Additionally, Rodrigue (2001) equation results
particularly inadequate to predict the existence of local
extrema on the typical terminal velocity curve.
Distilled water can still contain contaminants that
affect the rising motion of a single bubble. The
presence of such pollutants, especially surfactants,
changes the interface characteristics, rigidifying it
and reducing its shape oscillations during rise.
Duineveld (1995) and Sanada et al.
(2008)
determined experimentally the rise velocity of bubbles
moving through ultra-purified water with rectilinear
trajectories (20 < Re < 685). Fig. 2 shows
the comparison of those experimental data with
predictions of the proposed formulation and some
other from the literature, as the equation developed by
Tomiyama et al. (2002). The latter equation requires
the knowledge of the aspect ratio. For small equivalent
diameters it is still acceptable to relate this parameter

www.rmiq.org

Baz-Rodrguez et al./ Revista Mexicana de Ingeniera Qumica Vol. 11, No. 2 (2012) 269-278
to the Eo (Celata et al., 2007). Here, Sanada et al.
4 0

3 6

3 0

|V | ( c m / s )

2 0

1 8

(c m /s )

2 7

D u in e v e ld (1 9 9 5 )
S a n a d a e t a l. (2 0 0 8 )
E q . (1 ) w ith E q s . (6 ) a n d (8 )
F a n a n d T s u c h iy a (1 9 9 0 ) (%
J a m ia la h m a d i e t a l. (1 9 9 4 ) (%
T o m iy a m a e t a l. (2 0 0 2 ) (% E
R o d rig u e (2 0 0 1 ) (% E rr = 2 4

1 0

(% E rr
E rr = 2
E rr =
r r = 9 .2
.4 )

= 3 .2 )
1 .7 )
3 5 .1 )
)

A y b e rs a n d T a p u c u (1 9 6 9 )
E q . (1 ) w ith E q s . (6 ) a n d (8 )
9
0 .1

0
0 .0 4

0 .0 8

0 .1 2

d
e

0 .1 6

0 .2 0

(c m )

(2007) equation (Eo < 0.5) was used

1
1 + 6.5Eo1.925

0 .3

Fig. 2. Terminal velocity from experiments on

ultra-purified water (25 < Re < 685; 0.02 <
We < 3.36; 0.02 < Eo < 0.5) at 20-23 C.
Markers correspond to experimental data, and lines to
predictions from the equations.

E=

0 .2

(9)

In Fig. 2 one can observe that the relative average

error (percentage, % Err in the figures) is significantly
smaller for the formulation proposed in this work
(3.2 % absolute). Tomiyama et al. (2002) equation,
which is the next better fit, obtain a relative average
error of 9.2 %. This seems to confirm that viscous
effects are concentrated in a thin boundary layer for the
bubble size ranging from 0.04 to 0.18 cm of equivalent
diameter, according to the approximation from Eq.
(6).
As the bubble diameter increases, up to Eo about
0.5, the onset of path instability occurs (Duineveld,
1995; Sanada et al., 2008; Velduis et al., 2008),
and trajectories become zig-zagging or helical. The
corresponding bubble size range spans between the
local extrema of the typical terminal velocity behavior
(see Fig. 1). For pure water these extrema are between
de 0.15 0.65 cm, approximately (Mendelson,
1967; Loth, 2010), and according to the proposed
formulation between de 0.17 0.67 cm. These
local extrema are important because they indicate
significant transitions in the bubble rise mechanics.
The local maximum coincides with the path instability
onset, while the local minimum indicates the point
where inertial effects become dominant over surface
tension effects in Eq. (8).

0 .4

(c m )

Fig. 3. Comparison between predictions of Eq. (1)

and the single bubble velocity vector magnitude from
Aybers and Tapucu (1969).
For larger bubble size, the rise regime is dominated by
inertia, and bubbles adopt the characteristic spherical
cap shape, with the trajectory becoming approximately
rectilinear again.
For the terminal velocity between the local
extrema, three bubble diameters are possible. The
lowest diameter corresponds to a motion with
rectilinear path, the middle diameter is embedded in
the intermediate regime with oscillatory trajectory,
and the largest diameter corresponds to a preferably
inertial regime and spherical cap shape for the bubble.
In this sense, the description of the bubble motion
within the oscillating trajectory range is a very
complex task, and it remains a research topic in
recents numerical and experimental works (de Vries
et al., 2002; Mougin and Magnaudet, 2006; Shew and
Pinton, 2006a; Veldhuis et al., 2008).
Even if the predictions proposed in this work fit the
experimental data with pure water for many practical
applications (Re between 25 and 685) best than other
predictions, they overestimate slightly the terminal
velocity for larger bubbles. This occurs until the
inertial regime is reached, as shown in Fig. 1. This
overestimation might be caused by the presence of
contaminants in the water or by oscillatory trajectories.
Aybers and Tapucu (1969) measured both the
velocity magnitude and its axial component (terminal
velocity) for bubbles in the equivalent diameter range
between 0.1-1 cm following helical paths. Fig. 3
shows a comparison between the predictions of Eq.
(1) [with substitutions of Eqs. (6) and (8)] to the
obtained data for the velocity magnitude. The percent
relative error has an average of 5.5 %, as compared

www.rmiq.org

273

Baz-Rodrguez et al./ Revista Mexicana de Ingeniera Qumica Vol. 11, No. 2 (2012) 269-278

VT = |V| sin

(10)

where VT is the axial component, |V| is the velocity

magnitude and is the helix angle with respect to the
horizontal plane, which is given by
!
h
(11)
= tan1
DH
h corresponds to the vertical bubble displacement that
corresponds to a full oscillation cycle and DH is the
helix diameter. In terms of amplitude and frequency
of oscillation, one has
!
VT f
= tan1
(12)
2A
where f is the frequency and A is the amplitude.
Even if there are some investigations that explore
analytically the effect of the non-steadiness of the
flow in unconfined single rising bubbles, there are
still no adequate mathematical descriptions of the
bubble oscillatory motion and its perturbation effect
on the surrounding liquid. This is because the helix
angle determination in Eq. (10) had to be done from
empirical data.

274

Haberman and Morton (1953)

Peebles and Garber (1953)
Okazaki (1964)
Aybers and Tapucu (1969)
Blandn-Arrieta (1997)
Leifer et al. (2000)
Eq. (10) (% Err = 5.7)
Fan and Tsuchiya (1990) (% Err = 4.2)
Jamialahmadi et al. (1994) (% Err = 5.2)

3 2

(c m /s )

2 8

with experimental data by Aybers and Tapucu (1969).

This allows us to state that if the approximation fits
into Eq. (1) (with VT |V|), then its axial component
can be determined from the bubble trajectory. For the
case of small bubbles (de < 0.17) and large ones (de >
0.67) the velocity module corresponds approximately
to its axial component, or terminal, velocity. Yet
it is necessary to know the trajectory to determine
the velocity magnitude and its axial component for
oscillating motion.
Helical trajectories are much more common that
zig-zagging ones (Shew and Pinton, 2006a), and
correspond to the observed motion well above the
nozzle (Aybers and Tapucu, 1969). Zig-zagging
trajectories are associated to an unstable wake, so
these trajectories usually become helical at large
(Ellingsen and Risso, 2001; Shew and Pinton, 2006a).
As a result, it is advisable to assume that helical
trajectories characterize the whole range of oscillatory
motion, but with a stationary terminal velocity during
the bubble ascent (Shew and Pinton, 2006a; Shew et
al., 2006).
From the vectorial description of a helical path,
Aybers and Tapucu (1969) obtained the following
relationship between the axial velocity component and
its vectorial magnitude

2 4

2 0
0 .2

0 .3

0 .4

0 .5
e

0 .6

0 .7

(c m )

Fig. 4. Comparison between predictions of equations

and experimental data for bubble terminal velocity in
distilled water (560 < Re < 1500; 2.3 < We <
5.9; 0.5 < Eo < 5.9) a 20-28 C. Markers correspond
to experimental data, and lines to predictions from the
equations.
There are several investigations that studied the
oscillatory motion of single unconfined bubbles
(Aybers and Tapucu, 1969; Blandn-Arrieta, 1997;
Ford and Loth, 1998; Brucker, 1999; Ellingsen and
Risso, 2001; Wu and Gharib, 2002; Shew and Pinton,
2006a,b; Shew et al., 2006; Veldhuis et al., 2008). It
has been observed experimentally that the oscillations
amplitude A 0.0035 m and frequency f 5 Hz
(Loth, 2010). When these values are substituted into
Eq. (12) a correction for the terminal velocity is
obtained in Eq. (10). This correction depends on
the amplitude, frequency and experimental values for
the terminal velocity. The following equation was
consequently fitted
= 71.388 12.054de

(13)

Taking into account all the aforementioned facts, Fig.

4 plots the predictions of Eq. (10) considering the
helical paths followed by the bubbles. It should be
remarked that when de < 0.17 and de > 0.67, then
= 90 and the trajectory is rectilinear (sin = 1).
The corrected fit is similar to that of Fan and Tsuchiya
(1990) and Jamialahmadi et al. (1994) equations, all
with a relative average error of about 5 % (absolute).
Experimental data for other liquids than water are
scarce. Aqueous solutions and water are most relevant
for practical purposes. Even, Eq. (1), along with Eqs.
(6) and (8) can also be used to calculate the bubble
terminal velocity for pure liquids.

www.rmiq.org

Baz-Rodrguez et al./ Revista Mexicana de Ingeniera Qumica Vol. 11, No. 2 (2012) 269-278
Table 1. Critical equivalent diameter (dec ) for the onset of path instability of rising bubbles in pure liquids at
20-22 C.
Liquid

dec (cm)
[Eq. (1) with (6) and (8)]

dec (cm)
(Experimental data)

Reference

0.17
0.12
0.16
0.12
0.11
0.19

0.18
0.13
0.17
0.12
0.13
0.22

Duineveld (1995)
Hartunian and Sears (1957)
Hartunian and Sears (1957)
Hartunian and Sears (1957)
Zenit and Magnaudet (2008)
Zenit and Magnaudet (2008)

Pure water
Methanol
Ethanol
Benzene
Silicone oil (DMS-T00)
Silicone oil (DMS-T02)

8 0

6 0

R e la tiv e e r r o r (% )

Conclusions

A n ilin e
E th y l E th e r
A c e to n e
M e th a n o l
B e n z e n e
T o lu e n e
E th y l A c e ta te
Is o p ro p y l A lc o h o l
A c e tic A c id
M in e ra l O il
V a rso l
T u rp e n tin e
M e th a n o l

1 0 0

4 0

2 0

-2 0

0 .1

d
e

(c m )
1

Fig. 5. Relative error (%) for the predictions of the

proposed formulation compared to experimental data
for various pure liquids.

In Fig. 5 the relative errors for the formulation

predictions of this work with respect to experimental
data by Haberman and Morton (1953) and Peebles
and Garber (1953) are shown. Even if the error
values are not as small as for pure water, 72 % of
the points have a relative error less than 25 %. These
results are similar to those obtained by using Fan
and Tsuchiya (1990) and Jamialahmadi et al. (1994)
equations. With respect to the path instability critical
diameter prediction, the formulation of this work is
good for bubbles rising in pure liquids of low viscosity,
according to the comparison against experimental
results, shown in Table 1. For bubbles rising in viscous
liquids, it is possible that the oscillatory regime could
not be clearly observed, since it can be inhibited by
viscous effects. Anyway, more thorough studies for
pure liquids, other than water, are still necessary to
explain satisfactorily the motion of bubbles and their
corresponding regimes.

An equation was obtained for the terminal velocity

of bubbles rising in pure liquids, particularly valid
for the intermediate and inertial regimes of bubble
motion. The formulation uses a weighting of the force
balance obtained from boundary layer theory for nondistorted bubbles and an analytical equation coming
from a mechanical energy balance.
Predictions
obtained this way fit better the available experimental
data for bubbles rising through pure water (de >
0.04 cm), improving the results obtained by other
similar formulations that do not require bubble shape
geometrical parameters. Oscillatory motion is the
only situation that required the introduction of a
correction, assuming helical rising paths. In general,
the formulation proposed in this work showed to
be adequate to predict terminal velocity of rising
bubbles in pure liquids. Finally, the local extrema,
which bound the onset of the path instability (local
maximum) and the point where inertial effects begin
to dominate (local minimum), are well predicted by
the proposed formulation.

Acknowledgments
The authors respectfully acknowledge the Consejo
Nacional de Ciencia y Tecnologa (CONACyT Mexico) for their support through a Graduate
scholarship (SBR), as well as through Grant CB-2005C01-50379-Y.

www.rmiq.org

275

Baz-Rodrguez et al./ Revista Mexicana de Ingeniera Qumica Vol. 11, No. 2 (2012) 269-278

Nomenclature
amplitude of mean oscillation cycle for
helical bubble motion (m).
de
diameter of a spherical bubble of equivalent
volume (m).
dec
critical equivalent diameter for onset of
path instability in bubble motion (m).
DH helix mean diameter for spiraling path in
bubble motion (m).
E
aspect ratio between minor and mayor axes
of oblate ellipsoidal bubbles.
f
frequency of mean oscillation cycle for
helical bubble motion (Hz).
g
acceleration of gravity (m/s2 ).
h
vertical
bubble
displacement
that
corresponds to a mean oscillation cycle for
helical bubble motion (m).
|V|
magnitude of the bubble motion velocity
vector (m/s).
VT
terminal velocity of single bubbles (m/s).
VT 1 bubble rise velocity when viscous effects
are important (m/s).
VT 2 bubble rise velocity when surface tension
effects are important (m/s).
Greek letters
density difference between liquid and
gaseous phases (kg/m3 )
L dynamic viscosity of liquid phase (Pa s)
L density of liquid phase (kg/m3 ).

surface tension (N/m)

liquid phase velocity due to bubble

displacement (m/s).

helix mean angle for spiraling path in

bubble motion ( ).
Dimensionless numbers
Eo Eotvos number (db2 g/).
Re particle Reynolds number (L db VT /L ).
We Weber number (L db VT 2 /).

novel approach. Computers and Chemical

Engineering 25, 571-576.

References
Aybers, N. M. and Tapucu, A. (1969). The Motion of
gas bubbles rising through stagnant liquid. Heat
and Mass Transfer 2, 118-128.
Blandn-Arrieta, J. (1997).
Visualization and
Description of Bubble Dynamics (in Spanish).
Chemical Engineering MSc.
Thesis,
Universidad Autonoma Metropolitana Iztapalapa, Mexico.
Bozzano, G. and Dente, M. (2001). Shape and
terminal velocity of single bubble motion: a
276

Brucker, C. (1999). Structure and dynamics of the

wake of bubbles and its relevance for bubble
interaction. Physics of Fluids 11, 1781-1796.
Cai, Z., Bao, Y. and Gao, Z. (2010). Hydrodynamic
behavior of a single bubble rising in viscous
liquids.
Chinese Journal of Chemical
Engineering 18, 923-930.
Celata, G. P., DAnnibale, F., Di Marco, P., Memoli,
G. and Tomiyama, A. (2007). Measurements of
rising velocity of a small bubble in a stagnant
fluid in one and two component systems.
Experimental Thermal and Fluid Science 31,
609-623.
de Vries, A. W. G., Biesheuvel, A. and van
Wijngaarden, L. (2002). Notes on the path
and wake of a gas bubble rising in pure water.
International Journal of Multiphase Flow 28,
1823-1835.
Duangsuwan, W., Tuzun, U. and Sermon, P. A.
(2011). The dynamics of single air bubbles
and alcohol drops in sunflower oil at various
temperatures. AIChE Journal 57, 897-910.
Duineveld, P. C. (1995). The rise velocity and
shape of bubbles in pure water at high Reynolds
number. Journal of Fluid Mechanics 292, 325332.
Ellingsen, K. and Risso, F. (2001). On the rise of
an ellipsoidal bubble in water: Oscillatory paths
and liquid induced velocity. Journal of Fluid
Mechanics 440, 235-268.
Fan, L.-S. and Tsuchiya, K. (1990).
Bubble
Wake Dynamics in Liquids and Liquid Solid
Suspensions. Butterworth-Heinemann Series in
Chemical Engineering.
Ford, B. and Loth, E. (1998). Forces on ellipsoidal
bubbles in a turbulent shear layer. Physics of
Fluids 10, 178-188.
Haberman, W. L. and Morton, R. K. (1953). An
experimental investigation of the drag and shape
of air bubbles rising in various liquids. David
Taylor Model Basin Report 802, US Department
of the Navy, Washington, DC.

www.rmiq.org

Baz-Rodrguez et al./ Revista Mexicana de Ingeniera Qumica Vol. 11, No. 2 (2012) 269-278
Hadamard, J. S. (1911). Movement permanent
lent dune sph`ere liquide et visqueuse dans
un liquide visqueux.
Comptes Rendus
Hebdomadaires des Seances de V Academie
des Sciences (Paris) 152, 1735-1738.
Hartunian, R. A. and Sears, W. R. (1957).
On instability of small gas bubbles moving
uniformly in various liquids. Journal of Fluid
Mechanics 3, 27-47.

Mougin,
G. and Magnaudet,
J. (2006).
Wake-induced forces and torques on a
zigzagging/spiralling bubble. Journal of Fluid
Mechanics 567, 185-194.
Okawa, T., Tanaka, T., Kataoka, I. and Mori,
M. (2003).
Temperature effect on single
bubble rise characteristics in stagnant distilled
water. International Journal of Heat and Mass
Transfer 46, 903-913.

Jamialahmadi, M., Branch, C. and MullerSteinhagen, H. (1994). Terminal bubble rise

velocity in liquids. Chemical Engineering
Research and Design 72, 119-122.

Okazaki, S. (1964). The Velocity of Ascending air

bubbles in aqueous solutions of a surface active
substance and the life of the bubble on the same
solution. Bulletin of the Chemical Society of
Japan 37, 144-150.

Kulkarni, A. A. and Joshi, J. B. (2005). Bubble

formation and bubble rise velocity in gas-liquid
systems: A review. Industrial and Engineering
Chemistry Research 44, 5873-5931.

Peebles, F. N. and Garber, H. J. (1953). Studies on

the motion of gas bubbles in liquids. Chemical
Engineering Progress 49, 88-97.

Lehrer, I. H. (1976). A rational terminal velocity

equation for bubbles and drops at intermediate
and high Reynolds numbers.
Journal of
Chemical Engineering of Japan 9, 237-240.
Leifer, I., Patro, R. K. and Bowyer, P. (2000). A study
on the temperature variation of rise velocity for
large clean bubbles. Journal of Atmospheric and
Oceanic Technology 17, 1392-1402.
Levich, V. G. (1962).
Physicochemical
Hydrodynamics. Prentice Hall Englewood Clift.
Loth, E. (2010). Particles, Drops and Bubbles: Fluid
Dynamics and Numerical Methods. University
of Illinois at Urbana-Champaign & University
of Virginia.
Book Draft for Cambridge
University Press. Available on (May 7, 2010):
http://www.ae.illinois.edu/loth/CUP/Loth.htm.
Marks, C. H. (1973). Measurements of the terminal
velocity of bubbles rising in a chain. Journal of
Fluids Engineering 9, 17-22.
Mendelson, H. D. (1967). The Prediction of Bubble
Terminal Velocities from Wave Theory. AIChE
Journal 13, 250-253.
Moore, D. W. (1963). The boundary layer on
a spherical gas bubble.
Journal of Fluid
Mechanics 16, 161-176.
Moore, D. W. (1965). The velocity of distorted gas
bubbles in a liquid of small viscosity. Journal of
Fluid Mechanics 23, 749-766.

Rodrigue, D. (2001). Generalized correlation for

bubble motion. AIChE Journal 47, 39-44.
Ruzicka, M. C., Zahradnik, J., Drahos, J. and
Thomas, N. H. (2001).
Homogeneusheterogeneous regime transition in bubble
columns. Chemical Engineering Science 56,
4609-4626.
Shah, AND. T., Kelkar, B. G., Godbole, S. P.,
Deckwer, W-D., 1982. Design parameters
estimation for bubble column reactors. AIChE
Journal 28, 353-379.
Shew, W. L. and Pinton, J.-F. (2006a). Dynamical
model of bubble path instability. Physical
Review Letters 97(144508), 1-4.
Shew, W. L. and Pinton, J.-F. (2006b). Viscoelastic
effects on the dynamics of a rising bubble.
Journal of Statistical Mechanics: Theory and
Experiment P01009.
Shew, W. L., Poncet, S. and Pinton, J.-F. (2006).
Force measurements on rising bubbles. Journal
of Fluid Mechanics 569, 51-60.
Sanada, T., Sugihara, K., Shirota, M. and Watanabe,
M. (2008). Motion and drag of a single
bubble in super-purified water. Fluid Dynamics
Research 40, 534-545.
Talaia, M. A. R. (2007). Terminal velocity of a
bubble rise in a liquid column. Proceedings of
the World Academy of Science, Engineering and
Technology 22, 264-268.

www.rmiq.org

277

Baz-Rodrguez et al./ Revista Mexicana de Ingeniera Qumica Vol. 11, No. 2 (2012) 269-278
Tomiyama, A., Celata, G. P., Hosokawa, S. and
Yoshida, S. (2002).
Terminal velocity of
single bubbles in surface tension force dominant
regime. International Journal of Multiphase
Flow 28, 1497-1519.
Veldhuis, C., Biesheuvel, A. and van Wijngaarden,
L. (2008). Shape oscillations on bubbles rising
in clean and tap water. Physics of Fluids
20(040705), 1-12.
Wellek, R. M., Agrawal, A. K. and Skelland, A. H. P.
(1966). Shape of liquid drops moving in liquid
media. AIChE Journal 12, 854-862.

278

Wu, M. and Gharib, W. (2002). Experimental studies

on the shape and path of small air bubbles rising
in clean water. Physics of Fluids 14, L49-L52.
Zenit, R. and Magnaudet, J. (2008).
Path
instability of rising spheroidal air bubbles: A
shape controlled process. Physics of Fluids
20(061702), 1-4.
Zudin, AND. B. (1995). Calculation of the rise
velocity of large gas bubbles. Journal of
Engineering Physics and Thermophysics 68, 1015.

www.rmiq.org