Group Title: 7th International Conference on Multiphase Flow - ICMF 2010 Proceedings
Title: 9.4.2 - Investigation of mass transfer in poly-disperse gas-liquid systems by using a multi-variate population balance
Full Citation
Permanent Link:
 Material Information
Title: 9.4.2 - Investigation of mass transfer in poly-disperse gas-liquid systems by using a multi-variate population balance Computational Techniques for Multiphase Flows
Series Title: 7th International Conference on Multiphase Flow - ICMF 2010 Proceedings
Physical Description: Conference Papers
Creator: Buffo, A.
Marchisio, D.L.
Petitti, M.
Vanni, M.
Mancini, N.
Podenzani, F.
Publisher: International Conference on Multiphase Flow (ICMF)
Publication Date: June 4, 2010
Subject: population balance modelling
gas-liquid stirred tanks
quadrature-based moment methods
Abstract: The investigation of mass transfer in gas-liquid systems requires a formulation of a multivariate population balance equation, in order to consider the effect of a inhomogeneous distribution of bubble size, velocity and composition. A very efficient method to solve this equation is represented by the direct quadrature method of moments (DQMOM), where the evolution of the properties of the bubbles is evaluated by tracking several moments of the Number Density Function (NDF), with the advantage that the obtained equations can be easily implemented in commercial Computational Fluid Dynamics (CFD) codes. In this work, a detailed model was developed, tested on a simplified case study and coupled with a CFD code with the aim to provide a complete description of transport phenomena in a real system. The air-water stirred tank investigated by Laakkonen et al. (2006) was here simulated and results were compared with the experimental data.
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: VID00227
Source Institution: University of Florida
Holding Location: University of Florida
Rights Management: All rights reserved by the source institution and holding location.
Resource Identifier: 942-Buffo-ICMF2010.pdf

Full Text

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

Investigation of mass transfer in poly-disperse gas-liquid systems by using a
multi-variate population balance

A. Buffo* D.L.MarchisioT M.Petitti* M.VanniT N.Mancinitand F.Podenzanit

Dip. Scienza Mat. Ing. Chimica, Politecnico di Torino, Duca degli Abruzzi 24, 10129, Torino, Italy
t Division Refining & Marketing, Eni S.p.a., Via Felice Maritano 26, 20097, San Donato Milanese, Italy and
Keywords: population balance modelling, gas-liquid stirred tanks, quadrature-based moment methods


The investigation of mass transfer in gas-liquid systems requires a formulation of a multivariate population balance
equation, in order to consider the effect of a inhomogeneous distribution of bubble size, velocity and composition. A
very efficient method to solve this equation is represented by the direct quadrature method of moments (DQMOM),
where the evolution of the properties of the bubbles is evaluated by tracking several moments of the Number
Density Function (NDF), with the advantage that the obtained equations can be easily implemented in commercial
Computational Fluid Dynamics (CFD) codes. In this work, a detailed model was developed, tested on a simplified
case study and coupled with a CFD code with the aim to provide a complete description of transport phenomena in a
real system. The air-water stirred tank investigated by Laakkonen et al. (2006) was here simulated and results were
compared with the experimental data.


Agitated gas-liquid reactors are widely used in the
chemical and biochemical processes. In these systems,
gas dispersion is very important since it strongly influ-
ences gas-liquid mass transfer. Mechanical agitation en-
hances the homogeneity of the system improving the
overall performance of the reactor and flexibility of the
The design and scale-up of agitated gas-liquid tanks
is difficult, because chemical reactions usually have a
complex relationship with mass, momentum and heat
transfer; moreover, controlling phenomena often vary
with the reactor size. Traditionally, some semi-empirical
correlations based on vast experimentation were used
for the scale-up of reactors. These correlations are ex-
tremely useful for an initial analysis, but they consider
only volume-averaged properties calculated over the en-
tire system (i.e., global gas hold-up) and they can be
applied only to specific vessel geometries. However,
large industrial stirred tank reactors are often strongly
inhomogeneous: ,igilii k.ii spatial inhomogeneities ex-
ist even in laboratory scale tanks (Calderbank 1958)
and their effect on global properties must be consid-
ered. Computational Fluid Dynamics (CFD) has become
a popular tool for the analysis of chemical reactors, be-

cause it allows for the estimation of the fluid flow and
turbulent fields in almost any vessel size and geome-
try. At the moment, multiphase CFD predictions are un-
certain and require further development and validation
(Ranade 2002). CFD codes, together with the Popula-
tion Balance Models (PBM) and rigorous mass transfer
models, offer a great potential for the detailed analysis
of local mass transfer in agitated gas-liquid reactors.
As it is well known, gas-liquid mass transfer is a com-
mon rate-limiting step of reactor performance. Mass
transfer is often controlled by the liquid side through
the volumetric mass transfer coefficient (kLa). Corre-
lations for this parameter are popular, because the mea-
surement and modelling of the global mass transfer area
for unit volume is quite difficult (Calderbank 1958). Lo-
cal interfacial area for gas-liquid mass transfer depends
on the Bubble Size Distribution (BSD), which is known
to vary in agitated tanks depending on the operating con-
ditions and for the same operating condition from point
to point in the vessel. For example, bubbles near the
impeller break-up due to the high shear experienced,
whereas bubbles accumulate and coalesce in stagnant
zones. Moreover, every single bubble has a different
concentration of chemical species; this is important in
order to determine the driving force for the local mass
transfer flux. Population balances are a generalized ap-

proach for modelling the evolution of the local BSD, the
local distribution of concentration of the different chem-
ical species and the local bubble velocity distribution.
Of course to solve the PBM one must describe bubble
breakage and coalescence by using phenomenological
models. The final result is the formulation of a multi-
variate PBM capable of predicting both size and con-
centration distributions.
In the present work such an approach is formulated,
discussed and implemented in a commercial CFD code.
The approach is then used to model bubble coalescence
and breakage as well as mass-transfer. Model predic-
tions are validated by comparison through experimental
data from Laakkonen et al. (2006, 2007).

PBM for gas-liquid systems

A poly-disperse gas-liquid system can be successfully
described by using PBM. It is possible to consider a pop-
ulation of bubbles characterized by their size L, their
velocity Ub and their chemical composition 4b (inter-
nal coordinates), defining a Number Density Function
(NDF) so that the following quantity:

n(L, Ub, ?b; x, t) dL dUb d4b dx,

represents the expected number of bubbles, in the in-
finitesimal volume dx dx 1 dx2 dx3 around the phys-
ical point x (= (x, 2, X3) with size in the range be-
tween L and L + dL, velocity in between Ub and Ub +
dUb and composition in between 4kb and 4b + d4b.
The evolution of this NDF is governed by the Gener-
alized Population Balance Equation (GPBE) that reads
as follows:
On a a a
-t + a- (Ubin) + (Gn) + I (Abin) +
at 0dxi dL OUbi
+- f .) '9h(L, Ub, b;x,t), (1)

where UbV is the ith component of the bubble veloc-
ity, Ub (with i E 1, 2, 3), G is the rate of continuous
change of bubble size usually related to molecular pro-
cesses, where bubbles continuously grow (or shrink) be-
cause of the addition (or depletion) of single molecules
(e.g., evaporation or condensation processes), Abi is the
ith component of the acceleration, Ab, acting on a bub-
ble as a result of continuous forces, Ybj is the continu-
ous rate of change of bubble composition with respect
to the chemical component j (with j E 1,... q) as a re-
sult of chemical reactions and/or mass transfer between
phases. The term on the right-hand side is the discon-
tinuous event term, accounting in this case for instan-
taneous changes of size, velocity, and composition due
to bubble collisions, coalescence and breakage. In this

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

work it is assumed that G and 4bj depend only on mass
transfer of the different chemical components. More-
over, the acceleration acting on a bubble, Ab, is sim-
ply considered as a result of buoyancy, gravity and drag
force; any other inferfacial force is in this context ne-
glected, but could of course be added if necessary.
In the literature different methods are available to
solve the GPBE, but not all of them are able to treat
multivariate problems and not all of them result in
computationally efficient methods capable of describ-
ing real systems with complex geometries. For example,
both Classes Methods (CMs) and Monte Carlo Methods
(Vanni 2000; Zucca et al. 2007) are generally not used
for these problems, due to their computational requests.
It was demonstrated that methods based on quadrature
approximations, such as the quadrature method of mo-
ments (QMOM) (McGraw 1997; Marchisio et al. 2003),
and the direct quadrature method of moments (DQ-
MOM) (Fan et al. 2004; Marchisio et al. 2005), could
be amenable for solving PBM coupled with CFD codes;
QMOM has been recently used for the simulation of gas-
liquid systems (Petitti et al. 2010) however some limita-
tions where in that case detected. First of all the fact
that, all the bubbles are moving with same velocity, does
not allow for accounting for some segregation phenom-
ena that can be described only by tracking the distribu-
tion of bubble velocity. Secondly the introduction of ad-
dictional internal coordinates is tricky and not straight-
forward. For this specific application is interesting to
explore the possibilities offered by DQMOM. DQMOM
involves the direct solution of the transport equations for
weights and abscissas of the quadrature approximation.
Moreover, each node of the quadrature approximation
can be treated as a distinct gas phase, with its own ve-
locity, size and composition.
Applying the quadrature approximation, the NDF be-

n(L, Ub, kb; x,t)
= E wa(L L,)65(Ub Ub,a)6(b 0b,a),

where Wa is the number density (i.e., number of bub-
bles per unit of total volume) of the bubbles with size
equal to La, velocity equal to Ub,a and with composi-
tion vector equal to 4b,. Each delta function is char-
acterized by a different index a and it can be thought
of as a group of bubbles having the same size, veloc-
ity, and composition, each group representing one of N
nodes of the quadrature approximation. In order to re-
duce the number of equations, it is possible to consider
only two nodes (N = 2) with an acceptable lack of ac-
curacy as explained by Marchisio et al. (2005); more-

over, only two chemical components (resulting therefore
in only one chemical composition) are considered in this
work. Readers familiar with the method will recognize
that this is not a serious limitation since in DQMOM
internal coordinates can be easily added. For brevity in
what follows, space and time dependencies will be omit-
By substituting Eq. (2) into Eq. (1), and applying a
moment transformation, it is possible to derive the fol-
lowing transport equations:

a a
-(wi)+ ,(w1Ubi,l)
at ore,

al, (3)

at (w2) + a (W2Ubi,2)= a2,

a a
(wiLi) + a (W1Ubi,lL1) b +

SwiGi(L1, Ub,, 1b,1),

a a
at (w2L2) + x (2Ubi,22) = b+
+ w2G2(L2, Ub,2, Ob,2),

a a
a (Wlb,1) + a (W1Ubi,lb,1) = di

+ Wlbl(L1, Ub,1, b,1),

S(w24b,2) + 2 (2bi,2 b,2) = da2

+W2b2 (L2, Ub,2, Pb,2),

a a
a W (WlUbj,1) (W 1 Ubi,1 Ubj,l)
a~t Oxzi


+ lAbj(L, Ub,1, b,1), (9)

(w2Ubj,2) + (W2Ubi,2Ubj,2) = Cj2+
at Oxi
+W2Abj(L2, Ub,2, b,2), (10)

where the subscripts 1 and 2 refer to the two nodes of
the quadrature approximation.
Since only mass transfer of a single component is con-
sidered, Yb,i can be written as:

bi b(Li, Ubi, b,i) ( Hb,i), (11)

where kL is the mass transfer coefficient, 0c is the con-
centration of the chemical specie in the continuous liq-
uid phase and Ob,i is the concentration of the chemical

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

species in the gas bubble and H is the Henry constant.
Mass transfer affects also the size of the bubbles through
the bubble growth rate G. Through a simple mass bal-
ance on a single bubble it is possible to write that:

i = G(Li, Ubi, Yb,i) 2kLM- (0c Hb,i), (12)
where Mw is the molecular weight of the chemical
species considered. The acceleration acting on a single
bubble, Ab, can instead be written as follows:

ebAbi = bAb(Li, Ubi, Pb,i)
3 acCD,ic
4 L,

-0- + bS+

Ubi (U Ub),

where g is the gravity acceleration, yc and eb are the
densities of the continuous (liquid) and disperse (gas)
(5) phases, ac is the liquid volume fraction and (Uc Ubi)
is the slip velocity, namely the velocity difference be-
tween the continuous phase and the bubbles and CD,i
CD(Li) is the drag coefficient. There are several equa-
tions available in the literature to estimate CD,i; in this
work the Haberman and Morton model (reported in Wal-
lis 1969):
4 (Vc eb)Lg
CD,i 3 -U,2 (14)
3 PcUt2
(7) is used. This model is based on the terminal velocity of
a bubble of diameter Li in a liquid at rest, Ut, where is
evaluated through the Clift correlation (reported in Clift
et al. 1978):

( (2.14a
Ut = + 0.505 gLi)
V ecL

Li > 1.3mm

where a is the surface tension. Moreover an ad-hoc cor-
rection taking into account the effects of turbulence and
of the presence of other bubbles is also implemented in
this work (details are reported in Petitti et al. 2010). The
terminal velocity used to compute the drag force was
corrected with respect to the value obtained for an iso-
lated air bubble rising in water at rest, in order to take
into account for the effect of turbulence. Although this
correction should be applied and calculated on the local
turbulence intensity and gas hold-up, in this work a con-
servative approach was used, considering a constant ter-
minal velocity of 13 cm/s through out the entire stirred
The source terms of this set of equations (Eqs. 3-10)
are calculated by forcing the method to track some spe-
cific moments of the NDF. The choice of the moments
is completely arbitrary and in this work only pure mo-
ments are selected. This choice allows to decouple the

evolution of the different internal coordinates, as it will
become clearer below. It is examplificative to examine
how the source terms are calculated. Let us start with a1,
a*, b* and b*. These four quantities are calculated all to-
gether by solving the following linear system, written in
matrix form as:
Aac d, (16)
where the coefficient matrix stands for:
1 1 0 0

2 2 (17)
A -L2 2Li 2L2 (17)
-2L -2L2 3Li 3LJ

and the vector of unknowns a is defined by:

S= [al, a2, b", b ]T (18)


2kg M
b = bi wi ( Hb,i). (19)
The known right hand vector d can be written as follows:

d = [Sooo, S1oo, S2o, 30 oo]T (20)

where Skoo with k e 0, 1, 2, 3 corresponds to the source
term for the k* pure moment with respect to the size L.
By employing the quadrature approximation, this can be
written in the following closed form:

13 N N3
Skoo 2= W w (Li + L) /a(Li, L,)+
i=1 j=l
Lw, 5 a(Li, Lj)wj+
i=1 j=l

+ b(L .(kO)
+ f = (1 '0


where a(Li, Lj) is the coalescence kernel, b(Li) is the
breakage kernel and b() is the distributive daughter
function. For the moment these quantities are assumed
to depend only on bubble size. However, bubble compo-
sition dependencies can be easily accounted for.
Repeating the same procedure for d* and d*, it is pos-
sible to define:

A [1 1 ]2
A [21 212 '

a =[d*, d ]T,
d =[Solo, So20 + a al+ ]2a21]T,

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

d = di -w H b,i), (25)
and Solo with I E 1, 2 corresponds to the source term for
the th pure moment with respect to the composition y:

Solo = 2 S wi a(L L )w J L + L 3 +
i=1 j=1 1 3
wi a(Li, Lj)wj+
i=1 j 1

+ b(Lj)wjb t')
+i= 1


The expression of bW') differs from b ko), as it will be
explained below.
The same procedure is applied to obtain the source
terms for the momentum equations c,1 and c ,2. These
terms account for the change in momentum of the two
bubble classes due to discontinuous events such as bub-
ble collisions. For dilute systems, however, this effect
is neglegible with respect to drag and buoyancy forces,
therefore it will, in this contest, be neglected.
It is worth mentioning that this choice (k E
0, 1, 2, 3; 1 E 1, 2) allows us to track the first pure four
moments of the NDF with respect of the bubble size and
the moments of order one and two with respect to bubble
composition and velocity.
Before calculating Eqs. 21 and 26, it is necessary to
define the expressions for the distributive daughter func-
tions, coalescence and breakage kernels. In this work,
the following breakage kernel will be used (Laakkonen
et al. 2006):

b(L) 6.0 e/3 x

x erfc 0.04 2/3L5/3 + 0.01 Cd13L4/


where p is the viscosity of the continuous phase and e
is the rate of energy dissipation.
To model the frequency of bubbles coalescence it is pos-
sible to refer to kernel expression of Prince et al. (1990):

a(A, L) = 0.88 /3(A + L)2(A2/3 + L2/3)/2(A, L),
(22) (28)
where the coalescence efficiency can be written by using
the model of Coulaloglou et al. (1977):

(23)109 t .
(24) q(A,L) exp -610'A ) (29)

Concerning the distributive daughter function, it can be
realistic to assume binary breakage (as explained by An-
dersson et al. 2006) resulting in:

b(k,) 3240L (30)
(k + 9)(k + 12)(k + 15)
When 1 = 0, Eq. 30 becomes:

b(ko) 3240L (31)
(k + 9)(k + 12)(k + 15)'
and when k = 0, Eq. 30 can be written:

b') = 2y. (32)

The mass transfer coefficient kL is estimated by using
the Lamont and Scott model (Lamont et al. 1970):

k= 0.4 D5 o ,25 (33)

where D is molecular diffusivity of the chemical species
in the liquid phase and v, is the kinematics viscosity.
Now all the transport equations constituting DQMOM
are closed and through their solution it is possible to re-
construct the evolution of this poly-disperse gas-liquid
system. Coupling these equations to a commercial CFD
code, however, requires to rewrite them in terms of vol-
ume fractions:
7 3
ai =-Li wi, (34)

where ac is the volume fraction of the generic ith dis-
perse phase. Thus, the final governing equations be-
a a
at (alb) + -a (alebbi,l) 3

=e(b (2L- L al, (35)

a a
S(29eb) + (a2Pebbi,2)

a a
at (albL1) + a (a1ebLUbi,) )
at oxi

/(2 ( ,
Qb 7rL3 b

a a
a (02bL2) + (02PbL2Ubi,2)
at 9xi

2 La) ,

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

S(alebb,l) + 0, (aebbb,lUbi,1)

Qb (7 L1 + Yb,1 1b L a)) (39)

a a
(aQ2bb,2) + a (a2Pb-b,2Ubi,2)
Pb ( 3 2 2+ b,2 ( 2 Lb2 L- 22)) (40)
(6 2 2 (2 2 2 2 2L* 7 7

a a
a (cal bUbjI) + a (alpbUbi,lUbj,l)
a~t Oxzi

S-- + alebgj+
3 ocCDec
+ 1 IU- Ubj,1 (Ucj
4 L,

Ubj,1), (41)

a (a2PbUbj,2) + -a (2b bUbi,2Ubj,2)
= -a-2 + a2gbgi+
3 coCD c
+ 4 a D U Ubj,2 (Uj Ubj,2). (42)
4 L2
These equations must be solved together with the
equations for the continuous phase:
a a
at (ac&) + a- (ccUci) = 0 (43)

aa a xi
S(accUcj) + (accUciUcj) =
= --c- + acecgj+

+ ac(Pc + T,c) + )) +

2 OUck ,
3 xk
2 3 CD
4- ca k Uc Ub Icj bj,k)

S(Ocecc) + (cecQcOci)
at 6xi

,?c (kL6a (4c Hfb,j) + kL 60(c
L, L2


where PT,c is turbulence viscosity, modelled as:

2 L2 a2

Pb 2 L 3b

PT,c = 0.09 Pc

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

Table 1: Summary of gas and liquid properties used in
yc 998.2 kg m3
Pb 1.255 kg m3
fPc 1.0 10-3 Pas
Pb 1.78 -10-5 Pas
a 0.07 N m1

The values of turbulent kinetic energy k and rate of tur-
bulent energy dissipation c are obtained, in this work,
through k c Multiphase Mixture Model (Kataoka et al.

Test case, operative conditions and numerical

The gas-liquid stirred tank showed in Fig. 1 is investi-
gated in this work. The reactor is a 194 liter four-baffled
reactor, agitated by a six-blade Rushton turbine with cir-
cular metal porous sparger with diameter of 3.3 cm and
pores of 15 pm, located 10.5 cm below the impeller
and used for low gassing rates. This particular config-
uration was already investigated by Laakkonen et al.
(2006, 2007) and detailed experimental measurements
are available for model validation. The continuous liq-
uid phase is water and the disperse gaseous phase is air;
the properties used in the simulations are showed in Tab.
The simulations were carried out with the Multiple
Reference Frame approach and were performed by con-
sidering half of the geometry, since this is the smallest
portion of the reactor that satisfies the conditions of geo-
metrical symmetry and periodicity. 230,000 hexahedri-
cal cells were used to mesh the geometry, since both
the flow and the power number did not show, igMil Ik.iil
change with further grid refinement. The computational
grids were created by using GAMBIT, whereas simula-
tions were performed with the commercial code FLU-
ENT 12 and User-Defined Functions and Subroutines
were used to implement the DQMOM equations and the
drag force based on bubble terminal velocity. The trans-
port equations for DQMOM were implemented by defin-
ing additional scalars for the disperse gaseous phases.
CFD predictions were validated by comparison with ex-
perimental data. The predicted BSDs were compared
with the experimental ones in the measurement points
reported in Fig. 1 for a wide range of operating condi-
tions. The stirring rate ranged between 155 and 250 rpm
and the gas flow rate ranged from 0.018 to 0.093 volume
of gas per volume of reactor per minute (vvm), resulting
in global hold-up values up to 1.5 %. Validation through
comparison with transient mass transfer experiments (by

Figure 1: Rappresentation of stirred tank simulated and
location of sample points in which BSD data was de-

using oxygen and nitrogen) was also carried out and only
preliminary results are reported here.
As far as the boundary conditions are concerned, the
section of the sparger from which the gas is introduced
in the stirred tank was defined as a velocity inlet; the
gas volume fraction was evaluated as portion of the ge-
ometrical area available for gas flow and then the gas
velocity was determined by knowing the gas volumetric
flow rate. The liquid velocity components were set null.
The upper surface of the reactor was defined as a pres-
sure outlet in order to let the gas exit from the system
and, if there is a backflow, it is formed by liquid only.
Setting the boundary condition for the inlet bubbles it
is a non-trivial problem. The bubbles at the metal porous
sparger were modelled as a log-normal distribution as-
suming a standard deviation of 0.15, as usually observed
with metal porous spargers, and with a mean bubble size
at the inlet ( 1 . calculated with the correlation pro-
posed by Kazakis et al. (2008):

/. 7, 35 We- 1Re'1Fr1'8 dp171 (47)

where the dimensionless numbers of Reynols (Re), We-
ber (We), Froude (Fr) are defined as follows:

ecU2 d
We 9c U


and where U., is the gas superficial velocity based on
sparger area, d, is the sparger diameter and d, is the
mean pore size.

From this log-normal distribution the first mo-
ments with respect to bubble size were calculated and
then through the Product-Difference algorithm (Gordon
1968) weights and abscissas were also determined. With
this procedure the inlet volume fractions for the two bub-
ble classes (ai, a2) and their characteristic sizes (L1,
L2) were calculated. The inlet velocities were assumed
to be equal (Ub,i Ub,2 9.94 ms 1). The av-
erage inlet composition, written in terms of the inlet
oxygen concentrations in the two bubble classes (4b,1,
Yb,2), was instead fixed to be identical to that of the ex-
periments. For numerical and stability reasons also for
the composition a log-normal distribution with a stan-
dard deviation of 0.15 and a mean concentration of 8.56
mol/m3 was considered.
Some tests were also conducted on a homogeneous
system in order to verify the model and to detect all pos-
sible numerical problems. The set of equations was for
this simplified case solved in MATLAB.

Results and discussion

Firstly the tests carried out on the simple homogeneous
system will be presented, then the CFD simulations
and the comparison with experimental data will be dis-

a) PBM for homogeneous system. Several tests were
conducted using MATLAB for an hypothetical homoge-
neous system constituted of a simple computational cell.
In all of them mass transfer was not considered, and the
attention was focused on bubble coalescence and break-
up in order to evaluate the consistency of the multivariate
An example of these tests is reported in Fig. 2, which
shows the time evolution of the quadrature nodes and
weights in the case of N = 2 and in the case of pure
coalescence and pure breakage. As it is seen for pure
coalescence, the number of bubbles per unit volume,
wl and W2, decreases with time, while the characteristic
size of these two groups of bubbles, L1 and L2. When
pure breakage was simulated, the number of bubbles per
unit volume (the sum of wl and w ) increases with time,
while the characteristic sizes of these two groups of bub-
bles, L1 and L2, decrease.
In both cases the characteristic concentration of oxygen
in the two bubble classes, after a short period of time,
tends to the same value at about the average concentra-
tion. When such a situation takes place, the sizes of the
two classes are very different from the initial ones and
it is likely that the transfer of bubbles between the two
classes has been so large that the memory of the differ-
ent initial concentration has faded.
When simultaneous coalescence and breakage is simu-

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

lated with realistic kernels and volume-averaged turbu-
lent dissipation rate (c) (Fig. 3) the process is dominated
by coalescence. All the weights and abscissas of quadra-
ture approximation, wl, W2, L1, L2, Yb,i and Ob,2, as-
sume correct trends in time and physically consistent fi-
nal values. These tests showed the formal correctness
of the method that can therefore be implemented in the
CFD code.

Table 2: Comparison between experimental data and re-
sults from simulation.
Distribution of d32 (mm)
155 rpm and 0,018 vvm
Sample point Experimental Data Model
R2 2.37 3.10
R4 2.48 2.56
R8 2.29 2.57
R9 1.65 2.63
R12 3.31 3.09

Distribution of d32 (mm)
220 rpm and 0,041 vvm
Sample point Experimental Data Model
R2 2.56 2.66
R4 3.34 3.04
R8 2.57 2.47
R9 1.76 2.50
R12 3.81 3.20
Distribution of d32 (mm)
250 rpm and 0,052 vvm
Sample point Experimental Data Model
R2 2.74 2.45
R4 2.93 3.31
R8 2.17 2.55
R9 2.01 2.65
R12 3.18 3.57

b) CFD simulations on the real geometry. The re-
sults presented and discussed here mainly focus on the
prediction of BSD and their comparison with experi-
mental data from the literature (Laakkonen et al. 2006,
2007); the predictions of gas distribution profiles and
global hold-up for this configuration were the subject of
previous work (Petitti et al. 2010), and will be not dis-
cussed here.
An example of a typical simulation is reported in Fig.
4. As it is possible to see, the quadrature approxima-
tion, composed in this case of only two nodes (N = 2),
can be thought of as a description of the population of
bubbles in two classes. As the bubbles exit from the

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

10 x 3 15
?2 1 I
5 2 1 0.5

0 0 0 0
0 0.01 0.02 0.03 0 0.01 002 0.03 0 a01 002 0.03 0 0.01 0.02 0.03
1 t t t
0.02 0.03 1.5x10' 2 x107

0.02 I 1 2
j0.01 0,01 0.5 1 1

0 0 0 ----- 05
0 0.01 0.02 0.03 0 0.01 002 0.03 0 001 0.02 003 0 0.01 0.02 0.03
I t t t
10 10 10 10
95 95
S9 J9 > t 9 90
8.5 85
8 88 8
0 0.01 0.02 0.03 0.01 002 0.03 0 0.01 0.02 0.03 0 0.01 0.02 0.03
I t t t

Figure 2: Plots of weights and abscissas of NDF, calculated with DQMOM and kL 0, in the case of pure coalescence
(left) and pure breakage (right).



0 0.01 0.02 0.03 0.04 0.05

x 10


0 0.01 0.02 0.03 0.04 0.05


0.01 0.02 0.03 0.04 0.05

0 0.01 0.02 0.03 0.04 0.05


0 0.01 0.02 0.03 0.04 0.05

0.01 0.02 0.03 0.04 0.05

Figure 3: Plots of weights and abscissas of NDF, calculated with DQMOM and kL
coalescence and breakage.

0, in the case of simultaneous


7th International Conference on Multiphase Flow,

ICMF 2010, Tampa, FL, May 30 -June 4, 2010

i ,:. ,

S-- ,9.94et0D
9 54-.OD
9 24e+0

8.0 e+00
S' .... 7.45et00

6 194e+OD
96 1 e+00

5 "01 8 .5e O0
o, 0..0.537e-O0

SO O 417e00

u 4 t 7o .5e+00
2.31 6et00
2 09e+O
1 19e-00 1 393-00
4 95e01 795e-DO
5.96e-0O 4.6e 0
S00e00 017e+-D00
99 e+0D 9.604e+00
743e-00 8O 43e+O

8 40et00

528 e00, 1 337e+00
I .04e+00 I .0e+D00
4.57e -O I 546e-OO
4.32e+00 5 12e+00
3.8400 460e+DO
O+00 4.31e+00
S36-00 4.80e+DO
3.12-+00 3I7e+00

2 40eD 2 1.-et00
2.16e,00 7O63e.00
1680e00 2.00e300
1.44- 00 1.75e+00
1.20e+00 1.42e+DO
9.57.01 4.89e+00
477e-01 5574.01
2.37201 269e-01
2600e00 0 00-e000

0.018 vvm.

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






Operative condition
250 RPM 0.093 vvm

DfiuetLr (.am)

Figure 5: Comparison of Distribution of mean Sauter Diameter d32 (nun) and local Bubble Size Distribution.

sparger, they evolve due to coalescence, breakage and
mass transfer between liquid and gas, and their charac-
teristic volume fractions (ai, a2), sizes (L1, L2), veloci-
ties (Ub,1, Ub,2) and oxygen concentrations (0b,1, Yb,2)
change in tracking the evolution of the NDF as dictated
by the GPBE.
In Fig. 5, instead, the comparison between experi-
mental data and CFD predictions is reported for one op-
erative condition investigated. The comparison is here
limited to the BSD and the mean Sauter diameter. In
general, the model is able to describe qualitatively and
quantitatively the experimental trends also under other
operating conditions (see Tab. 2). The experimental
data show a clear trend for the Sauter diameter (d32): in
all the operative conditions it is possible to observe that
dR12 > dR4 > dR8 and dR2 > dRg, where the subscript
indicate the sampling points. Moreover, in general, d32
increases with the gas flow rate: in this case, in fact, the
bubbles collide and coalesce more frequently, whereas
the breakage frequency decreases. The values detected
in zone R2 and R8 are very similar, but for high flow
rates the d32 measured in R2 is higher than in R8; in
other words, for high flow rates, coalescence becomes
more important near the upper surface, where the bub-
bles tend to accumulate.
The comparison with CFD results show that the agree-
ment with the experimental data increases as the gas
flow rate and the rotational speed of the impeller are in-
creased. At 155 RPM, all the predicted values are larger
and the trend for the experiments is not clearly observ-

At 220 RPM, there is good agreement with the experi-
mental data; the trend for the Sauter diameter is clearly
observable: dR12 > dR4 > dR8 and dR2 > dRg, with
the increase of the distance between dR12 and dR9. For
this rotational speed and gas flow rate it is clear that the
bubbles with the smallest size stay on the impeller axis,
as the lower recirculation zone entraps the bubbles with
the largest size.
At 250 RPM, the general qualitative agreement is very
good. For 0.052 vvm, the plane near the impeller is not
well described: in fact slightly dR9 > dR2. Whereas for
0.093 vvm, the experimental trend is well represented
by the simulations, as confirmed by the progressive in-
crease of coalescence near the upper surface.

Let us now discuss the results concerning mass trans-
fer. In Fig. 6 the contour plot of the local volumetric
mass transfer coefficient, kL a, is reported over the whole
geometry. This variable is strongly inhomogeneous, in
fact the larger values are on the impeller plane and near
the buffles, where breakage caused by turbulence signif-
icantly increases the number of bubbles. In these zones
there is the largest interfacial area per unit of volume a
and the higher mass transfer coefficient kL. On the con-
trary, in the bottom part of reactor, where no bubbles are
present, the mass transfer coefficient nearly vanishes. It
is important to remind that most of these effects could
not be appreciated without using a population balance
model. Moreover with this approach not only the size
and the velocity distributions are correctly described, but

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






Figure 6: Spatial distribution of the local mass transfer
coefficient kLa (m s 1) for 250 rpm 0.093 vvm.

* Experimental Data [220 RPM 0,041 vm]
* Experimental Data [250 RPM 0,052 vm]
- CFD prediction [220 RPM 0,041 wm]
- CFD prediction [250 RPM 0,052 wm]

Time (s)


Figure 7: Comparision between CDF predictions and
experimental data of oxygen concentration YL (mol
m 3) considering a monodisperse bubbles distribution
kLa (m s 1) for different operative conditions.

the composition distribution is also accounted for. In
fact, by solving Eqs. 39 and 40, the method is capable
of modelling the fact that smaller bubbles are charac-
terized by higher mass transfer coefficients than bigger
Since the code still suffers of some numerical issues
that jeopardize its stability, in Fig. 7 only some prelim-
inar results of two mass transfer tests for a monodis-
perse bubble distribution are reported. These experi-
ments were carried out fluxing into the stirred tank ni-
trogen and then switching to oxygen. The time evolution
of species concentration in the liquid was tracked for a
single sample point (R4) (Laakkonen et al. 2006, 2007).
As it is seen there is indeed room for improvements.
Further steps of this work will focus on the solution of
all the numerical problems, related to the CFD imple-

In this work, a multivariate population balance for poly-
disperse gas-liquid systems was formulated. DQMOM
has been used to solve the PBM for a homogeneous sys-
tem to consider only the effect of coalescence and break-
age on the bubbles. Tests conducted on this simple ge-
ometry show the consistency of the model in the case of
pure coalescence or breakage as well as combined coa-
lescence and breakage problem.
Then the equations were rewritten in order to be cou-
pled with a commercial CFD code to simulate the evolu-
tion of the NDF for a real geometry, for which exper-
imental data concerning local BSD and mass transfer
were available. The comparison between CFD predic-
tions and experimental data for the real system indicates
generally good agreement, even using only two nodes
for the quadrature approximation. A first mass transfer
analysis was successful carried out but some numerical
issues did not allow a complete model validation. The
next steps of our work will focus on these issues and on
the final model validation.


Andersson R., Andersson B., On the breakup of fluid
particles in turbulent flows, AIChE Journal, Vol. 52, pp.
2020-2030, 2006

Calderbank P.H., Physical rate processes in industrial
fermentation, Part 1: The interfacial area in gas-liquid
contacting with mechanical agitation, Trans. Inst. Chem.
Engrs., Vol. 36, pp. 443-463, 1958

Clift R., Grace J.R., Weber M.E., Bubbles, Drops and
Particle, Academic Press, London, 1978

Coulaloglou C.A., Tavlarides L.L., Description of in-
teraction process in agitated liquid-liquid dispersion,
Chemical Engineering Science, Vol. 32, pp. 1289-1297,

Fan R., Marchisio D.L., Fox R.O., Application of the di-
rect quadrature method of moments to polydisperse gas-
solid fluidized beds, Powder Technology, Vol. 139, pp.
7-20, 2004

Gordon R.G., Error bounds in equilibrium statistical me-
chanics, Journal of Mathematical Physics, Vol. 9, pp.
655-667, 1968

Kataoka I., Serizawa A., Basic equations of turbulence
in gas-liquid two-phase flow, International Journal of
Multiphase Flow, Vol. 15, pp. 843-855, 1989

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

Kazakis N.A., Mouza A.A., Paras S.V., Experimental
study of bubble formation at metal porous spargers: ef-
fect of liquid properties and sparger characteristics on
the initial bubble size distribution, Chemical Engineer-
ing Journal, Vol. 137, pp. 265-281, 2008

Laakkonen M., Alopaeus V., Aittamaa J., Validation of
bubble breakage, coalescence and mass transfer mod-
els for gas-liquid dispersion in agitated vessel, Chemical
Engineering Science, Vol. 61, pp. 218-228, 2006

Laakkonen M., Moilanen P., Alopaeus V., Aittamaa J.,
Modelling local bubble size distributions in agitated ves-
sels, Chemical Engineering Science, Vol. 62, pp. 721-
740, 2007

Lamont J.C., Scott D.S., An eddy cell model of mass
transfer into surface of a turbulent liquid, AIChE Jour-
nal, Vol 16, pp. 513-519, 1970

Marchisio D.L., Virgil R.D., Fox R.O., Quadrature
method of moments for aggregation-breakage process,
Journal of Colloidal and Interface Science, Vol. 258, pp.
322-334, 2003

Marchisio D.L., Fox R.O., Solution of population bal-
ance equations using the direct quadrature method of
moments, Journal of Aerosol Science, Vol. 36, pp. 43-
73, 2005

McGraw R., Description of aerosol dynamics by the
quadrature methods of moments, Aerosol Science &
Technology, Vol. 27, pp. 255-265, 1997

Petitti M., Nasuti A., Marchisio D.L., Vanni M., Baldi
G., Mancini N., Podenzani F., Bubble size distribution
in stirred gas-liquid reactors with QMOM augmented by
a new correction algorithm, AIChE Journal, Vol. 56, pp.
36-53, 2010

Prince M.J., Blanch H.W., Bubble coalescence and
break-up in air-sparged bubble columns, AIChE Journal,
Vol. 36, pp. 1485-1499, 1990

Ranade V.V., Computational flow modelling for chem-
ical reactor engineering, Academic Press, San Diego,
USA, 2002

Vanni M., Approximate population balance equations
for aggregation-breakage processes, Journal of Colloid
and Interface Science, Vol. 221, pp. 143-160, 2000

Wallis G.B., One Dimensional Two-phase Flow, Mc
Graw Hill, New York, 1969

Zucca A., Marchisio D.L., Vanni M., Barresi A.A., Val-
idation of bivariate DQMOM for nanoparticle processes
simulation, AIChE Journal, Vol. 53, pp. 918-931, 2007

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

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