7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
Numerical study of the mass/heat transfer from deformed bubbles rising in stationary
liquid
Bernardo Figueroa, Dominique Legendre
University de Toulouse, IMFT / INPTUPS France
All6e du Prof. Camille Soula, Toulouse, 31400, France
bfigueroae@yahoo.com, legendre@imft.fr
Keywords: mass transfer, bubble deformation, aspect ratio, boundary fitted grid
Abstract
Mass (or heat) transfer was studied for the case of an ellipsoidal bubble rising through a stationary viscous liquid. The numerical
code JADIM that solves the NavierStokes equations, coupled with Diffusionadvection of a passive scalar was used to
characterize the effect of the bubble deformation on the interfacial transfer. Different types of boundary fitted numerical grids
were tested in order to obtain reliable results for the mass/heat transfer (Adoua et al, 2009). The theoretical predictions based on
potential flow theory (Lochiel and Calderbank, 1964, Favelukis, 2005) present some differences with real flows due to the
vorticity boundary layer wake and structure and a correction depending on the Reynolds number is thus necessary. Existing
correlations for the Sherwood/Nusselt number have been tested using the bubble equivalent radius as characteristic length
(Takemura and Yabe, 1998). These correlations seem valid for intermediate to large Reynolds numbers of order 102 to103.
Introduction
Mass transfer in gasliquid flows is a very complex problem
of recognized industrial and academic relevance. The present
investigation was carried over within the context of a
multidisciplinary project aimed to develop a suit of
computational tools to simulate and design aeration basins in
water treatment plants1. The specific task of determining the
effect of bubble deformation on the transfer was investigated
numerically, in order to provide correlations used as closure
model for large scale Eulerian simulations. Previous
investigations (Moore, 1963, Duineveld, 1995, Riboux et al.
2010) have shown that the shape of millimetricsized bubbles
rising in water can be reasonably approximated by oblate
spheroids, of equivalent radius
r= (ab2)1 (1)
where a and b are the minor and major semiaxes of the
ellipsoid of aspect ratio
X= a/b (2)
used to describe the bubble deformation. This simple
geometry allowed for the theoretical study of mass/heat
transfer to the surrounding medium (Lochiel & Calderbank
1964), in the limit of large bubble Reynolds numbers
2r Uoo
Req = (3)
V
where U, is the bubble terminal velocity and V is the
kinematic viscosity of the liquid phase. Mass transfer at the
bubble surface is normalized using the Sherwood number
see https://o2star.cemagref.fr/leprojetxxx for more
information
(analog to the Nusselt number in heat transfer)
fDVc'dA
Sh, A (4)
S DVc'A
L
where D is the diffusion coefficient, c' is the concentration of
some gaseous species (for example oxygen) and L is a
characteristic length. When L=2req the subscript 'eq' will be
used. In what follows, whenever the major semiaxis b is
used as characteristic length, no subscript will be added.
When the advectiondiffusion equation is brought to
dimensionless form, the Peclet and Schmidt numbers appear
as relevant parameters:
Sc = ; PeL = ReLS (5)
D
The subscript L will be specified depending on the
corresponding characteristic length taken for each particular
case. The equivalence between the Reynolds and Sherwood
numbers, based on b and req are then
Req = Reb 13
(6)
Sheq =Shb 1/3
For gas bubbles, mass or heat transfer is driven by the
external flow hydrodynamics, the flow field inside the
bubble being negligible due to the large viscosity difference
between the liquid and the gas (at least two orders of
magnitude).
Some theoretical results are available in the literature for
both the steadystate and transient regimes for spheres
(Boussinesq, 1905, Ruckenstein, 1959) and ellipsoidal
inclusions (Lochiel & Calderbanks, 1964, Favelukis, 2005).
If the Reynolds number is large enough, then it can be
assumed that the external flow in the liquid phase is potential,
with the exception of a thin boundary layer at the bubble
surface. For the particular case of spheres, in the limit
Re oo (Boussinesq, 1905):
Shs (Re > o)= Pe/2 (7)
Some experimental investigations can also be found on the
subject for oblate and prolate spheroids (see for example
Skelland & Cornish, 1963) where empirical correlations are
also obtained in terms of Reynolds, Peclet and Schmidt
numbers. A broad bibliographical review for momentum and
Mass/heat transfer with and without phase transition can be
found in Michaelides i 'i, I'i).
In a more recent investigation, Takemura and Yabe (1998)
calculated the oxygen concentration of millimetric (almost
spherical) gas bubbles by means of precise partial pressure
and (timevarying) diameter measurements. They proposed a
semiempirical expression to calculate the Sherwood number,
which will be henceforth referred to as TY and used later for
comparison purposes.
TY = 1 2 3 (2.5 + Pe12)
TY = 3 (1+0.09Re2/)3/41
(8)
A description of the methods and geometrical details are
given in the next section.
Implementation and validation
The code used is the 'JADIM' code developed at the Institut
of Fluid Mechanics of Toulouse. It can solve the unsteady
NavierStokes equations in terms of velocitypressure
variables, for any orthogonal 3D curvilinear coordinate
system. The discretization method is finite volumes and the
precision is second order in time and space using
RungeKutta/CrankNicolson schemes (see Magnaudet,
Rivero & Fabre, 1995, Calmet & Magnaudet, 1997 and
Legendre & Magnaudet, 1998). This code has been used in
several studies concerning hydrodynamics and transfer
problems for bubbles and particles.
Boundary fitted grids are produced by choosing a coordinate
transformation that forces the transformed curvilinear
coordinates to match the interface geometry (while keeping
the other components constant). In the particular case of
ellipsoids, the curvilinear 2D ellipsoidal system described in
Adoua et al. (2009) was used to generate the geometry. This
kind of configuration is advantageous because it allows for a
very good mesh refinement near the interface, which is
necessary if thin thermal and hydrodynamic boundary layers
are to be described. Figure 1 depicts one example of this
kind of grid. The total domain size is 80 times the major
semiaxis b, in order to minimize confinement effect due to
the finite domain size. In the upper right of the figure a
zoomin of the bubble is shown. The liquid velocity at
infinity is represented by an arrow showing the direction of
the undisturbed liquid flow. The whole xy plane is turned
around the x axis, so as to produce an axisymmetric domain,
with an oblate spheroidal bubble at the origin. Two types of
grids were tested, the second one was similar to that of
Legendre & Magnaudet (1998), where the coordinates were
generated from the equipotentials and streamlines of the flow
around ellipsoidal cylinders. The results obtained with these
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
two grids are very closed and we decided to the first kind of
grids (ellipsoidal coordinates).
I '
120
80
50
50
Figure 1: an example of 'ellipsoidal' grid. The 'radial'
coordinate far from the bubble are almost spherical.
The range of parameters of the simulations was chosen such
as the flow field remains axisymmetric: if a given threshold
of Reynolds number and eccentricity is exceeded, there is a
transition (path instability), as described by Mougin &
Magnaudet ii" ,). Under these conditions vortex shedding
is produced by vorticity accumulation in the bubble equator,
and an oscillating lift force causes the bubble trajectory to
destabilize. Figure 2 shows a map of the Re vs X plane.
Each symbol can represent several simulations since
different Sc were considered. The curve to the left represents
the transition to an axisymmetric wake, as calculated by
Mougin & Magnaudet i,2i"I ). The curve to the right
represents the transition to a 3D wake or path instability. The
dotted lines are approximations of the original curves, of the
form:
X (Re) = 0.182(ln Re)3 + 0.391(ln Re)2
2.5341nRe +6.826 (9)
2 (Re) = 0.029(ln Re)3 + 0.685(ln Re)2
5.202 In Re + 14.970
The arrows point to illustrations of the corresponding type of
flow regime. The filled ellipses are located at the
approximate locus of an air bubble in water whose diameter
is indicated by the scale at the top of the figure (deq in mm).
Note that the crosses correspond to bubbles with
axisymmetric wakes, while the asterisks have 3D wakes. In
what follows, results will be presented for bubbles with
axisymmetric flow fields unless stated otherwise.
1.4 1.8 2.4 2.8
U 0
10 00 0 0 D
0 0
0 0 0
1 15 2 2.5 3
Figure 2: transitions to axisymmetric and 3D wakes.
Several validation tests were carried over with different
values of the parameters for spherical bubbles. Different grid
parameters, Reynolds and Schmidt numbers and finally
different aspect ratios were used to compare the
hydrodynamics and mass transfer with previous
investigations. The tests showed that convergence is more
difficult as the Peclet number increases, since the
hydrodynamic and thermal boundary layers get thinner and a
more refined grid must be used to describe the physics
correctly. Once a suitable set of grid parameters was chosen,
bubbles were simulated for different values of aspect ratio,
Reynolds and Schmidt numbers. The preliminary results are
presented in Figure 3, in the ShPe plane. The agreement
with TY is very satisfactory.
104
0 Sc=1
a Sc= 100
y Sc 300
O< Sc = 500
10' TYSc1 /
4, TY Sc 100
TY Sc = 300
TY Sc = 500
 2/1 11Pe ,/2
hydody' i an .ra onar aesgt hne n
10 ..... .. ........ ........ ........ .....
10 11 10' 101 102 103 1e0 10s 106
Pe
Figure 3: Shs vs Pe. Comparison with previous
investigations for spherical bubbles.
In the next section, the results of similar choices for Re and Sc
are presented using different aspect ratios varying from 1.15
to 3. The differences are discussed in the light of the details
obtained from the local transfer and velocities at the surface.
Numerical Results
Let us now consider Figure 4, which presents the local
Sherwood number at the bubble surface (normalized with
respect to the asymptotic value of Boussinesq, 1905) as a
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
function of the meridional angle: to the upper left, the local
transfer for an aspect ratio of 1.2 is presented for different Pe,
as well as the theoretical result of Favelukis i' 11i). The
agreement between theory and simulations is good in this
case: indeed, the differences can be explained in terms of the
velocity boundary layer (Moore, 1963, 1965). As the bubble
deformation increases, the structure of the boundary layer is
more complicated, and this velocity difference is more
important due to the effect of vorticity. As a consequence the
real transfer is smaller than the theoretical predictions, as can
be observed from the calculations for more deformed bubbles
(see Figure 4 (b), (c) and (d)). Note in (c) that the
recirculation zone has a non negligible effect (see the arrow
near the horizontal axis of (c), where the onset of the
recirculation is recognized), yet this transfer decreases with
time until the recirculation zone attains a very small transfer
rate, as hypothesized by Lochiel & Calderbank (1964). The
details of the time it takes for this to happen is currently being
investigated.
(a) X=12
(c) =2
25
1.5
I \
0.5
onsetof '
0 w rcirhon '
0 50 100 150 200
(b})I1.5
S.Pe=1 DODsim
1.5 P e=IDOO sim
 Pe=9DOOD sim
1 theory
1 l F=aulis2005
0
0 50 100 150 200
(d) 2.5
2.5
1.5
05 5 1 ..1
0 50 100 150 200
Figure 4. Local Sherwood number for different aspect ratio.
Figure 5 shows the ratio between the Sherwood number of
an ellipsoid and that of a sphere \!i 14I calculated using the
theory of Lochiel and Calderbank (1964), that depends on the
aspect ratio only (dashed line). The horizontal axis
corresponds to the aspect ratio, the different markers
representing different values of Sc. The agreement is good
for large Reynolds number, as expected from the theoretical
(potential flow) assumptions. For large deformations, some
discrepancies are evident. Apparently the transfer increases
for strongly deformed bubbles. This is an artificial effect
coming from the use of the major axis as characteristic length,
as will be discussed later. As already stated from the results
in Figure 4, one should expect the real transfer to be smaller
than the theoretical asymptotic resulting from potential flow
theory.
Sh/Sh forSc==I
sf.r
Figure 5: Sh based on b for different Sc and X
Note that for Re>300 the transfer becomes almost
independent of the Reynolds number and the curves are
almost superposed.
Figure 6 shows the simulated Sherwood number as a
function of the aspect ratio for different Sc. The markers
correspond to different Reynolds numbers (based on b). The
dotted line corresponds to an empirical equation obtained
from the simulation data:
Shb 1/(0.0129Re+2.278)1/3 for Re Re*
Sh, (10)
{ 1/(0.0129Re*+2.278)1/3 for Re > Re*
where Re*=350. This constant value represents the
behaviour of the curves at large Reynolds, where the transfer
becomes almost independent of this parameter. If Re* is not
introduced, the error starts to become important for values of
the Reynolds number above 1000, approximately.
Sc=10 S 100
4 1 2 3 1 5 2 2 
  
11Q'[t& __0 ^Sell 11&PGEBE
(t& ^  ^i ___ ___ _o O'( e   <   
1 15 2 25 3 1 15 2 25 3
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
normalizes the volume), i.e. if we compare two bubbles with
same Reeq and Peq based on the equivalent diameter and not
on their major axis, the transfer is seen to decrease sensibly
with the deformation. Figure 7 shows the Sherwood number
based on the equivalent diameter 1/, for Sc=500 and aspect
ratios x of 1, 2 and 3. The solid line corresponds to
Boussinesq (1905) solution. The circles represent spheres
(X=l), while the squares are ellipsoids of aspect ratio X=3 for
the same Reynolds number based on b. The corresponding
equivalent Peclet and Reynolds numbers do not correspond
anymore (see the diagonal arrow from the circle to the
square). TY is presented as a dashed line.
S=500
Figure 7: Sherwood number 1, / as a function of Peeq
If the formula TY is considered, but instead of using the
original arguments, one introduces the equivalent Sherwood
and Reynolds numbers in TY:
TYe =TY(Req,Scq) (11)
The equivalent Sherwood formula (11) can make the points
to collapse in a single curve, as shown in Figure 7. Note that
two points with the same Reynolds based on b and different
aspect ratio do not have the same Reeq. This 'collapse' was
also observed for the other Schmidt numbers considered.
I ___ sc=lo a R =10
x R =10
gg Ot R 500
/ o < "
Ii 2 2C
1 15 2 2S
M
Sc=100
101 2 1
1 15 2 26 3
S Se=500
10 2 2
1 *I ; T X ,L
11 2 25 3 10 '2 25
Figure 6: Sherwood number as a function of aspect ratio for
different Sc.
However, this parameterization is 'tricky', in the sense that
the transfer is normalized with respect to the bubble area and
not the volume of the ellipsoids. Indeed, if one takes a
Sherwood number based on the equivalent diameter (which
Figure 8. Equivalent Sherwood number for different Reeq.
Figure 8 shows the same kind of results as Figure 6, but in
this case the equivalent Sherwood number is plotted for
different values of equivalent Reynolds, represented as
dotted lines with different markers. TYeq corresponds to the
bullets. It is clear that it fits satisfactorily the simulation data
for large equivalent Reynolds, as expected since TY has the
same structure as the equations obtained by Lochiel &
Calderbanks (1964) for the limit of large Reeq. Both formulas
Sh/Sh for Sc==10
s her
.h
12
0.5 1 15 2 2.5 3 3.5
SWSh lorSc300
''1 ..x
s, ,
1 X X
.r .
0.5 1 1.5 2 2.5 3 3.5
ShLhl
{ Re=10
X Re.100
X Re=30
0 Re300
A RMOOO
05 I 1.5 2 25 3 3.5
Sh/Sh h a lor Sc= 00
1 ]
0.5 I 1.5 2 2.5 3 3.5
ShSh h forSc50
svhe"
sr=Vo0
(10 and 11) give reasonable errors (less than 10%) if
compared with the simulated data, relation (11) being more
suited for calculations in the context of water treatment
facilities, since the volume of gas injected into the tanks gives
the input power necessary to run the process and
consequently the efficiency. However eq. (10) behaves
slightly better, particularly if small Reynolds are also
considered.
Discussion
It is interesting to consider the existence of an optimum
bubble diameter in terms of transfer efficiency (for a given
volume of gas injected). Small bubbles have larger residence
times, but smaller transfer coefficients, while larger bubbles
have better transfer with smaller residence times. In order to
explore this question let us consider first, for simplicity, the
case of spherical bubbles:
To first order,
Sh ~ Pe12 r/2U2 (12)
The amount of matter per unit time is
mr = nDA. 2.7, I, (13)
where n is the number of bubbles. On the other hand, the
total transfer depends on the resident time, i.e. on the
terminal velocity. The latter is attained when buoyancy and
drag balance each other so
F, r3 ~ FD ~ U r => U ~ r2 (14)
where FB and FD are buoyancy and drag forces. From (13)
and (11):
Sh r3/2 (15)
The residence time can be estimated from h/ U If we
define the efficiency /7 as the amount of gas transferred
divided by the amount of gas injected, then it can be
obtained by multiplying eq. (12) by the residence time and
dividing it by the amount of gas injected, which is
proportional to r ; so one has to conclude that the
efficiency 7 scales as
h Sh h 1 1
M = m 2~ (16)
UQ r U r /2 1/216)
so there is no maximum nor minimum: the smaller the
bubbles, the better efficiency. For very small bubbles, the
Reynolds number tends to zero and Sh tends to a constant
value, giving 77 r This problem has been addressed in
the past by Motarjemi and Jameson (1978), with similar
conclusions. The authors calculated a depletion distance,
defined as the distance travelled by the bubble in order to
deliver 95% of its oxygen content. The effect of the
variations in oxygen concentration with depth was taken
into account, and the depletion distance h95% was presented
as a function of the depth of release. It was concluded that
usual bubble sizes (larger than 2mm) do not allow for
efficient oxygen transfer since they leave the reactor before
delivering more than 40% of their oxygen. Most of the
actual facilities (with depths between 1 and 5m) would have
to produce bubble equivalent diameters in the range 300
1000itm to attain reasonable levels of efficiency.
A characteristic time can be obtained to describe the
exponential decay of the bubble concentration with time.
First let us not that the variation of the concentration inside
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
the bubble evolves as (the mass flow rate divided by the
bubble volume):
dco 3 ShD
dc (c c) (17)
dt 2 r,,
and the evolution in time of the oxygen concentration in the
bubble can be obtained in virtue of Henry's law:
c, = Hcco (18)
where He is Henry's constant dimensionlesss), s ,
co and c are the concentration in the liquid at the interface,
in the gas and in the liquid far from the bubble, respectively.
As a consequence, the characteristic time to transfer all the
oxygen is
2 r2
T e (19)
S 3 ShDHc
An effective tank height can be calculated using
h95% = 3zrgU
If bubbles leave the tank before 'giving up' their oxygen
(below h95) the transfer is not optimized. Note that it is
assumed that C remains constant throughout the process.
For a water tank and bubbles of reqlmm this gives about
5m. Even if this calculation is rather simplistic, this value is
close to the one calculated by Motarjemi and Jameson[28]
(1978). Discrepancies appear for larger bubble diameters
because in our calculation, the variations of the interface
concentrations with depth are not considered. This
characteristic time T is anyway a good indicator of the
order of magnitude for both the depletion time and height,
as shown in Figure 9, and may be used in dimensional
analysis to estimate the overall efficiency.
 iari.UIH,
..^.. J s"n 1",
x 1
 ^.
1D& 10,
h
Figure 9: depletion height vs. equivalent diameter:
comparison with more elaborated implementation of
Motarjemi and Jameson
Conclusions and perspectives
The effect of the deformation on the mass/heat transfer was
studied numerically, for different Schmidt and Reynolds
numbers that correspond roughly to the one observed in
millimetric bubbles rising in water (following rectilinear
paths). Some existing correlations were compared to the
simulations, and it was found that TY fitted the numerical
results satisfactorily. If the parameters are normalized with
respect to the equivalent radius, the curves collapse to a
single curve that could help to describe the effect of
deformation on the transfer for practical purposes. A brief
discussion on the overall efficiency of the transfer was also
presented, where it is evidenced that practical facilities
produced 'too large' bubbles. A simple correlation for the
depletion time (and consequently the effective depletion
length) was given, as a function of the bottom depth (of
bubble production). It would be interesting to calculate the
effect of surrounding bubbles (induced agitation) on the
concentration, as well as the effect of foreaft symmetry,
which is lost for large Reynolds numbers. Some of our results
were not presented here because it was out of the scope of
this report, however, it should be noted that we succeeded to
produce foreaft asymmetric grids, and it was observed that
the effect of this shape distortion can be treated as in Zenit &
Magnaudet i,21" 1), with a judicious renormalization. The
hydrodynamics, as well as the transfer are affected by the
path stability for an ellipsoidal bubble, which depends
strongly on the aspect ratio and Reynolds number (Mougin &
Magnaudet, 2007). It is common to find air bubbles in water
with equivalent diameters between 2 and 4mm and aspect
ratios between 1.2 and 3mm (see filled symbols in Figure 2).
Additional 3D simulations are necessary to study this regime.
Acknowledgements
This work was granted the support of the National Research
Agency of France (Agence Nationale de la Recherche),
Project 02STAR, Reference: ANR07ECOT00701. The
authors would also like to thank their colleagues from the
CNRS Federation FERMaT, the IMFT, the LISBP and the
Cemagref d'Antony for making this collaboration possible.
References
Adoua, S.R., Legendre, D. & Magnaudet, J. Reversal of the
lift force on an oblate bubble in a weakly viscous linear shear
flow. J. Fluid Mech., 628, 2341 ,21 ',
Boussinesq, J. Title. J. Math. Pure Appl., 6, 285 (1905).
Calmet, I. & Magnaudet, J. Largeeddy simulation of high
Schmidt number mass transfer in a turbulent channel flow,
Phys. Fluids, 9, 118 (1997).
Favelukis, M. & Hung Ly, C. Unsteady mass transfer around
spheroidal drops in potential flow, Chem. Eng. Sci., 60,
70117021 i'"a 1).
Legendre, D. Quelques aspects des forces hydrodynamiques
et des transferts de chaleur sur une bulle sphY'erique, PhD
Dissertation, INPT, Toulouse, France (1996).
Legendre, D. & Magnaudet, J. The lift force on a spherical
bubble in a viscous linear shear flow, J. Fluid Mech, 368,
81126(1998).
Lochiel, A.C. & Calderbank, PH. Mass transfer in the
continuous phase around axisymmetric bodies of revolution.,
Chem. Eng. Sci., 19, 471484 (1964).
Magnaudet, J., Rivero, M. and Fabre, J. Accelerated flows
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
past a rigid sphere or a spherical bubble. Part 1. Steady
straining flow. J. Fluid Mech., 284, 97135 (1995).
Michaelides, E.E. Hydrodynamic Force and Heat/Mass
Transfer From Particles, Bubbles, and Drops The Freeman
Sholar Lecture, J. Fluids Eng., 125, 209238 (2003).
Moore, D. W. The boundary layer on a spherical gas bubble, J.
Fluid Mech., 16, 161176 (1963).
Moore, D. W. The velocity of rise of distorted gas
bubble in a liquid of small viscosity, J. Fluid Mech., 23,
749 (1965).
Motarjemi, M. and Jameson, M. Mass transfer from very
small bubbles the optimum bubble size for aeration. Chem.
Eng. Science, 33, 14151423 (1978).
Mougin, G. & Magnaudet, J. Wake instability of a fixed
spheroidal bubble, J. Fluid Mech., 572, 311337 ,2" ).
G Riboux, F. Risso & D. legendre 2010 Experimental
characterization of the agitation generated by bubbles rising
at high Reynolds number, J. Fluid Mech. 643, 509539.
Ruckenstein, E. On heat transfer between vapour bubbles in
motion and the boiling liquid from which they are generated,
Chem. Eng. Sci, 10, 2230 (1959).
Skelland, A.H.P & Cornish, A.R.H. Mass transfer from
spheroids to an air stream, A.I.Ch.E. Journal, 9/1, 7376
(1963).
Takemura F & Yabe, A. Gas dissolution process of spherical
rising bubbles, Chem. Eng. Sci., 53/15, 26912699 (1998).
