Chemical Engineering Journal xxx (2012) xxx–xxx
Contents lists available at SciVerse ScienceDirect
Chemical Engineering Journal
journal homepage: www.elsevier.com/locate/cej
A hybrid CFD framework for ﬂuidized bed ozonation reactors coupling
interface tracking and discrete particle methods
Rodrigo J.G. Lopes ⇑, Rosa M. Quinta-Ferreira
Centro de Investigação em Engenharia dos Processos Químicos e Produtos da Floresta (CIEPQPF), GERSE – Group on Environmental, Reaction and Separation Engineering,
Department of Chemical Engineering, University of Coimbra, Rua Sílvio Lima, Polo II – Pinhal de Marrocos, 3030-790 Coimbra, Portugal
h i g h l i g h t s
g r a p h i c a l a b s t r a c t
" An IT-DP framework was
Interstitial ﬂow snapshots of hydrodynamic ﬂow patterns resulting from IT-DP simulations for bubble
plume structures and bubble development (uG = 0.1 cm/s, [O3] = 20 g/Nm3, T = 20 °C, and P = 1 atm).
"
"
"
"
successfully developed for the
catalytic ozonation of phenol-like
pollutants.
The effect of physicochemical
properties was investigated in terms
of bubble diameter and detachment
time proﬁles.
The hybrid CFD model handled
agreeably the inﬂuence of ozone
velocity on multiphase ﬂow
patterns.
The interstitial ﬂow mappings of
ozone concentration unveiled a
heterogeneous degree of transport
phenomena.
The mineralization efﬁciency was
intimately correlated with bubble
plume hydrodynamic structures.
a r t i c l e
i n f o
Article history:
Available online xxxx
Keywords:
Interface tracking
Discrete particle method
Fluidized bed reactor
Catalytic ozonation
Hydrodynamics
Reactive ﬂow
.
a b s t r a c t
In this work, a ﬂuidized bed reactor was investigated theoretically and experimentally for the catalytic
ozonation of phenol-like pollutants. First, an interface tracking sub-model dealing with the individual
motion of ozone bubbles in the continuous phase, and a discrete particle sub-model accounting for the
trajectories of the solid particles were embedded accordingly into a hybrid CFD framework. Second,
the effect of physicochemical properties including the surface tension and ﬂuid viscosity was thoroughly
evaluated at different hydrodynamic ﬂow regimes. Afterwards, the inﬂuence of ozone velocity was quantiﬁed both on the gas and liquid superﬁcial velocities and on the detoxiﬁcation efﬁciency of liquid pollutants. Here, the interface tracking-discrete particle hybrid model was found to slightly overestimate
radial bubble velocities being almost negligible in the center region of the ﬂuidized bed reactor. As the
surface tension increased, the mineralization efﬁciency was considerably lower and became higher when
low-viscosity conditions were used for the ozonation of organic pollutants. Finally, the morphological
features of interstitial ﬂow maps at different hydrodynamic and reactive catalytic ozonation conditions
highlighted the occurrence of bubble plume-like structures which affected the overall decontamination
efﬁciency of phenol-like pollutants within the gas–liquid–solid ozonation reactor.
Ó 2012 Elsevier B.V. All rights reserved.
1. Introduction
⇑ Corresponding author. Tel.: +351 239798723; fax: +351 239798703.
E-mail address: [email protected] (R.J.G. Lopes).
Ubiquitous applications of three-phase ﬂuidized bed reactors
vary from chemical and petrochemical industries to biochemical
1385-8947/$ - see front matter Ó 2012 Elsevier B.V. All rights reserved.
http://dx.doi.org/10.1016/j.cej.2012.06.143
Please cite this article in press as: R.J.G. Lopes, R.M. Quinta-Ferreira, A hybrid CFD framework for ﬂuidized bed ozonation reactors coupling interface tracking and discrete particle methods, Chem. Eng. J. (2012), http://dx.doi.org/10.1016/j.cej.2012.06.143
2
R.J.G. Lopes, R.M. Quinta-Ferreira / Chemical Engineering Journal xxx (2012) xxx–xxx
Nomenclature
A
C
D
D
Eö
F
G
k’
I
m
_
m
p
R
r
Re
Si
Sc
Sm
Sh
interfacial area, m2
specie mass concentration, ppm
diameter, m
diffusivity, m2 s1
2
Eötvös number, Eö = (qL qb)gdb =r, dimensionless
force vector, N
gravity acceleration, ms2
reaction rate constant, min1
unit tensor, dimensionless
mass, kg
mass transfer, kg s1
pressure, Nm2
radius, m
reaction rate
Reynolds number, Re = qL|uB uL|db/lL, dimensionless
source term in the species balance equation
Schmidt number, Sc = lL/(qLD), dimensionless
volume averaged momentum transfer due to interphase
forces
Sherwood number, Sh = kLdb/D, dimensionless
and environmentally-based technologies. In the ﬁrst group, they
have been often used for oil upgrading processes, Fischer–Tropsch
synthesis, coal liquefaction and gasiﬁcation, whereas recently they
have been employed for bio-oxidation processes for waste gases
and wastewaters [1–3]. The multiphase ﬂow regimes typically
encountered in ﬂuidized beds have been reviewed elsewhere
addressing the hydrodynamic operation often carried out by a
cocurrent gas and liquid upﬂow system with liquid as the continuous phase [4]. Conversely to gas–solid systems, investigation of
gas–liquid-solid ﬂuidized beds has focused primarily on the nonreactive ﬂow operation both from an experimental and theoretical
viewpoint. The majority of these studies were concerned with the
intrinsic disadvantages of ﬂuidized beds such as solids backmixing,
attrition of particles, and erosion of surfaces aiming to achieve ultimately an uniform and regular ﬂow distribution. However, the reliable design and scale-up have to account for not only the
hydrodynamic description of ﬂow regimes, but also the accurate
replication of reactive ﬂow patterns.
Here, multiphase computational ﬂuid dynamics (CFDs) provides
realistic two- and three-dimensional representations of hydrodynamic and reaction parameters by means of three different frameworks for the investigation of ﬂuidized beds: Eulerian–Eulerian (E–
E) method, Eulerian–Lagrangian (E–L) method and direct numerical simulation (DNS) method [5–7]. While the Euler–Euler models
represent the two phases as interpenetrating continua [8,9], the
Euler–Lagrange models describe the liquid phase as a continuum
and track each individual bubble as it rises due to buoyancy forces
[10,11]. Euler-Euler and Euler–Lagrange models intrinsically
engender a coarse-grid representation of the bubbly ﬂow as the
spatial resolution down to the Kolmogorov scale became computationally intractable. Both theoretical frameworks allow merely
capturing the meso- and macroscale features of the multiphase
ﬂow ﬁeld so the numerical simulations of reactive and mixing phenomena are representatively-validated for characteristic reaction
times in the order of the meso- and macromixing temporal scales
[12,13].
In case of fast chemical reactions, sub-grid closure laws have to
be developed for enhancing the accuracy of the coarse-scale models as long as the temporal scale simulations demand high-resolution computations of microscale ﬂow. Aiming to predict the
interface position and motion, three different techniques have
t
uL
ub
V
time, s
liquid velocity vector, ms1
bubble velocity vector, ms1
volume, m3
Greek letters
a
volume fraction, dimensionless
l
viscosity, kg m1 s1
q
density, kg m3
r
interfacial tension, Nm1
s
stress tensor, Nm2
Subscripts
b
bubble
eff
effective
i
ith species
P
pressure
s
superﬁcial
T
turbulent
been employed including the moving-grid, the grid-free and the
ﬁxed-grid method. Speciﬁcally, the moving interface problem with
ﬁxed and regular grids is numerically solved by the front tracking
and front capturing methods [5,14]. The interface position is
tracked explicitly in the front tracking approach by the advection
of the Lagrangian markers on a ﬁxed/regular grid in contrast with
the implicit representation of the moving interface in the Eulerian
treatment by a scalar-indicator function deﬁned on a ﬁxed/regular
mesh point. The advection equation of the scalar-indicator function is solved for capturing the movement of the interface, which
is reconstructed by piecewise segments at every time step. Accordingly, the surface tension force is accounted for as a source term
using the continuum surface force method [15] so three subsequent methods have been proposed including the VOF method
[16], the marker density function, and the level-set method [17].
Speciﬁcally, the front tracking approach enables the direct and
accurate computation of the surface tension force by avoiding
the calculation of the interface curvature. Indeed, the Lagrangian
representation of the interface circumvents the interface reconstruction from the local distribution of the volumetric phase fractions, and overcomes the artiﬁcial merging of interfaces, as faced
by lattice Boltzmann and VOF frameworks.
In the realm of environmentally-based applications of ﬂuidized
bed reactors, the present case-study addresses the long-standing
interest in catalytic ozonation of high-strength liquid pollutants.
The modeling approach for this gas–liquid–solid reaction system
is twofold: to employ the front tracking sub-model for simulating
the individual motion of ozone bubbles in the continuous phase
consisting of polluted water and suspended catalyst particles, and
to use the discrete particle sub-model to account for the trajectories
of the solid particles by including the external forces and non-ideal
particle–particle and particle–wall collision events. This methodology allows investigating thoroughly the bubble swarm behavior
due to the interrelating dependency of translation and rotation of
ozone bubbles, and tackles the characteristic gravity, pressure, drag,
lift and virtual mass forces acting on the catalytic particles.
2. Previous work
A CFD–VOF–DPM method has been used by Li et al. [18] to perform numerical simulations of gas–liquid–solid ﬂuidization
Please cite this article in press as: R.J.G. Lopes, R.M. Quinta-Ferreira, A hybrid CFD framework for ﬂuidized bed ozonation reactors coupling interface tracking and discrete particle methods, Chem. Eng. J. (2012), http://dx.doi.org/10.1016/j.cej.2012.06.143
R.J.G. Lopes, R.M. Quinta-Ferreira / Chemical Engineering Journal xxx (2012) xxx–xxx
systems. They have combined the discrete particle method and a
volume-tracking approach to investigate the bubble wake behavior. Yang et al. [19] have presented a theoretical study of particle–ﬂuid systems by analysing two-ﬂuid CFD models. Critical
assertions were pointed out for these methods to gather the
meso-scale heterogeneity at the scale of computational cells,
which lead to inaccurate calculation of the interaction force between particles and ﬂuids. Regarding the homogeneous and bubbling ﬂuidization regimes, Renzo and Maio [20] have reported
DEM–CFD simulations for evaluating the hydrodynamic stability
of gas and liquid ﬂuidized beds. The appearance of bubbles in the
ﬂuidized bed behavior was demonstrated to occur at velocities in
quantitative agreement with the theory of ﬂuidized bed stability.
The transition of homogeneous ﬂuidized beds to bubbling ﬂow
regime has been investigated by means of Euler–Euler model
for monodisperse suspensions of solid particles. Based on this
approach, Mazzei and Lettieri [21] have evaluated both the
steady-state expansion proﬁles of liquid-ﬂuidized systems and
the stability of homogeneous gas-ﬂuidized suspensions.
Later, a three dimensional transient model was developed by
Panneerselvam et al. [22] for simulating the local hydrodynamics
of a gas–liquid-solid three-phase ﬂuidized bed reactor. A good
agreement with the experimental data has been claimed regarding
the ﬂow ﬁeld predicted by the CFD framework. Aiming to examine
the residence time distribution of a three-phase inverse ﬂuidized
bed, Montastruc et al. [23] compared CFD simulations with phenomenological semi-empirical models. They found that the increase of the gas ﬂow rate led to higher mixing intensity of the
gas phase, while the liquid phase performed closer to disperse plug
ﬂow. Sivaguru et al. [24] have presented several hydrodynamic
studies on three-phase ﬂuidized beds. The liquid and solid ﬂow
was represented by the mixture model, the air was treated by
means of discrete phase method, and they have found that computed pressure drop agreed well with the experimental data at different ﬂuid density conditions. A three dimensional CFD model of
the riser section of a CFB have been reported by Behjat [25] for
investigating the catalyst particle hydrodynamic and heat transfer.
An Eulerian model was used to model both gas and catalyst particle phases, whereas a Lagrangian approach was employed to simulate the ﬂow ﬁeld and evaporating liquid droplet characteristics.
For liquid–solid ﬂuidization systems, Huang [26] has evaluated
the effect of drag correlation on the multiphase ﬂow hydrodynamics. A drag correlation was proposed according to the CFD simulations also accomplished for the added mass force variable. Roghair
et al. [27] have reported a comprehensive investigation on the drag
force of bubbles in bubble swarms at intermediate and high Reynolds numbers. The authors have thoroughly developed an advanced front-tracking model for investigating bubble swarms
based on direct numerical simulations. Both operating and physiochemical variables were examined to gather a novel drag correlation needed for larger-scale models.
To the best of our knowledge, scarce CFD simulations have been
used so far to investigate a gas–liquid–solid system under reactive
ﬂow conditions. In this regard, the catalytic ozonation of liquid pollutants is comprehensively evaluated both from a theoretical and
experimental viewpoint. The remainder of this paper is organized
as follows. First, the interface tracking (IT) approach is embedded
with a discrete particle (DP) method to simulate the ﬂuidized
bed reactor. The multiphase model is presented with the constitutive equations and the boundary conditions. The experimental procedure is described after the simulation setup. Second, the
succeeding section addresses the ﬁrst results regarding the effect
of physicochemical properties including the surface tension and
ﬂuid viscosity. Then, the inﬂuence of ozone velocity is evaluated
both on the gas and liquid superﬁcial velocities and on the detoxiﬁcation efﬁciency of liquid pollutants. Finally, the morphological
3
features of interstitial ﬂow maps are discussed for validation activities at different hydrodynamic and reactive catalytic ozonation
conditions.
3. Mathematical model
The interface tracking-discrete particle framework encompasses two distinct models to simulate the multiphase ozonation
reactor. The interface tracking approach is used to compute the
motion of ozone bubbles in the polluted water. The continuous
phase is then characterized by the phenol-like compounds solubilized in the aqueous phase and by the suspended catalyst particles.
The discrete particle method numerically treats the motion of solid
particles by dealing with the overall action of external forces and
non-ideal particle–wall and particle–particle collisions.
3.1. Interface tracking model
Under unsteady conditions, the constitutive multiphase ﬂow
equations for Newtonian ﬂuids are expressed by
@ aL
þ r aL u ¼ 0
@t
@a u
q L þ r aL uu ¼ aL rp þ qg þ ðr aL l½ru þ ruT Þ
@t
þ Fr FLS
ð1Þ
ð2Þ
The source terms Fr and FLS represent the surface tension force
and two-way coupling due to the presence of the suspended catalyst particles, while aL expresses the volume fraction of aqueous
phase. The ﬂuid viscosity l and density q are calculated from the
local distribution of the phase indicator that is solved by the Poisson equation expressed in Eq. (3) as described elsewhere [28]. The
vector quantity G includes the local information of the spatial distribution of the gas–liquid–solid interface
r2 F ¼ r G
ð3Þ
The local density q and viscosity l is usually calculated according to Eqs (4) and (5) by considering the densities of the liquid and
gas phase.
q ¼ F qL þ ð1 FÞqG
l ¼ F lL þ ð1 FÞlG
ð4Þ
ð5Þ
However, the kinematic viscosities of liquid and gas phase have
been calculated by using a different method developed elsewhere
[29]
q
q
q
¼ F L þ ð1 FÞ G
l
lL
lG
ð6Þ
In the momentum balance equation, the surface tension force is
modeled by including the source term Fr . The direct computation
of the curvature of gas–liquid–solid interface is circumvented by
using the preliminary estimation of the net surface tension force
Fc;r acting on a speciﬁc surface node, c
Fc;r ¼
I
rðt nÞdc
ð7Þ
Eq. (7) requires for each separate edge e of the element c, to be
mapped to the Eulerian grid for obtaining the volumetric surface
tension force needed for the momentum equation. Typically, the
control volumes for the momentum balance comprise more than
one surface node so the individual node contributions are ultimately summed as described by
Please cite this article in press as: R.J.G. Lopes, R.M. Quinta-Ferreira, A hybrid CFD framework for ﬂuidized bed ozonation reactors coupling interface tracking and discrete particle methods, Chem. Eng. J. (2012), http://dx.doi.org/10.1016/j.cej.2012.06.143
4
Fr ðxÞ ¼
R.J.G. Lopes, R.M. Quinta-Ferreira / Chemical Engineering Journal xxx (2012) xxx–xxx
1
DxDy Dz
PP
c
q Dðx xc;k Þrðtc;k nc;k Þ
Pc;kP
c
k qc;k Dðx xc;k Þ
k
ð8Þ
The Newton’s second law has been used to deal with the individual motion of the suspended catalyst particles by accounting
for the pressure, gravity, lift, drag, and added mass forces
d
ðmp v Þ ¼ Fp þ Fg þ FL þ FD þ FVM
dt
¼ V p rp þ mp g C L V p qL ðv uÞ ðr uÞ
V p C D;eff
ðu v Þ
ð1 aL Þ
D
ðqL V p C VM ðv uÞ þ qL V p C VM ðv uÞ ruÞ
Dt
þ
ð9Þ
The force balance expressed in Eq. (9) is derived from the literature as to perform numerical simulations by means of a discrete
bubble model for dispersed gas–liquid two-phase ﬂow [30]. As
long as the size of the catalyst particles is considerably smaller
than the spatial resolution of the computational grid needed to
properly replicate the bubble dynamics, a supplementary relationship for the effective drag coefﬁcient (CD,eff) is required to close the
momentum balance equation for the solid phase. Here, the model
originally developed by Hoomans et al. [31] has been used to deal
with the drag between the aqueous phase and the catalyst particles. The effective drag coefﬁcient is expressed by taking into account two distinctive approaches: the Wen and Yu equation in
the dilute regime (aL > 0.8) and the Ergun equation in the dense regime (aL < 0.8) as shown respectively in Eqs. (10) and (11).
aL < 0:8 C D;eff ¼ 150
3.4. Catalytic ozonation kinetics
The catalytic ozonation kinetics of phenol-like compounds have
been derived considering a three-step mechanism with the compounds lumped into three groups: easier degraded pollutants (A),
intermediates with difﬁcult degradation (B) and the desired products (C). The reaction rates are given by Eqs. (14) and (15)
0
aL ð1 aL Þ
dp
qL ju v ja2:65
L
ð1 aL Þ2
lL
aL
2
dp
þ 1:75ð1 aL Þ
ð10Þ
qL
dp
ju v j
ð11Þ
CVM = 0.5 and CL = 0.5 were assumed during all IT-DP simulations, and the drag coefﬁcient CD was computed following wellknown correlations for drag coefﬁcients of spherical particles.
3.3. Species continuity equations
l
@
ðaL qL C L;i Þ þ r aL qL~
v L C L;i eff ;L rC L;i
@t
Sci
_
_
¼ mb;L C b;i mL;b C L;i þ aL Si
0
k3 C TOC B
ð14Þ
0
k2 C TOC A
ð15Þ
After numerical integration, the normalized total organic carbon (TOC) concentration is expressed as shown in Eq. (16)
0
0
0
0
0
0
C TOC
k2
k k
¼
ek3 t þ 0 1 0 3 0 eðk1 þk2 Þt
C TOC 0 k01 þ k02 k03
k1 þ k2 k3
ð16Þ
This model encompasses two different routes for reactivity of
0
oxidizable compounds: a direct conversion to end-products ðk1 Þ
and a ﬁnal oxidation preceded by formation of intermediates
0
0
ðk2 ; k3 Þ, and has been calibrated with the experimental conversion
data obtained in the range of inlet ozone gas concentration used in
our case-study. Table 1 summarizes the kinetic parameters of the
0
apparent reaction rate dependent on ozone k ¼ f ðC O3 Þ] for each
one of the three reaction steps involved in the application of ozone
inlet gas concentrations in the range 20–88 gO3/Nm3.
4. Numerical solution
4.1. Gas–liquid–solid ﬂow ﬁeld
A state-of-the-art ﬁnite volume technique has been employed
to solve the Navier–Stokes equations on a staggered rectangular
three-dimensional computational mesh as described by Lopes
and Quinta-Ferreira [32]. An implicit treatment of the pressure gradient coupling the two-step projection–correction method and an
explicit treatment of the convection and diffusion terms have been
used to numerically compute the gas–liquid–solid hydrodynamics
under reactive ozonation ﬂow conditions. The spatial discretization
of the convective terms was performed using a second-order ﬂux
delimited Barton-scheme, and a standard second-order central
ﬁnite difference scheme was employed for the diffusion terms.
Additionally, the incomplete Cholesky conjugate gradient algorithm was used to solve the pressure Poisson equation.
In what regards the density calculation method, a phase marker
function F for the respective spatial distribution has been implemented for the gas and solid phases. This phase indicator was
computed from the triangulated interface by solving a Poissonequation as described elsewhere [28]
r2 F ¼ r G ¼ r The species continuity equation is expressed in Eq (12), where
the concentration of chemical species i in the liquid phase is accounted by CL,i and Si is the source term describing the production
or consumption of species i due to homogeneous chemical
reaction.
0
r TOC A ¼ dC TOC A =dt ¼ ðk1 þ k2 ÞC TOC A
r TOC B ¼ dC TOC B =dt ¼
3.2. Discrete particle model
3
4
ð13Þ
i¼1
This approach enables an enhanced distribution of the surface
tension force by implementing commonly the dimensionless distribution function. In order to attain stable computations, Eq. (8) is further normalized by the cell volume for achieving the force density
function that is accounted for in the constitutive momentum equation. Apart from conferring additional and improved characteristics
for more prominent surface tension phenomena, this methodology
promotes a smoothening effect needed for the numerical simulation
of high-density ratio physicochemical systems.
aL > 0:8 C D;eff ¼ C D
ns
X
Ci ¼ 1
X
Dðx xc Þnc Dsc
ð17Þ
c
In Eq. (17), the summation is accomplished over all surface
nodes c representing the gas–liquid–solid interface, and expressing
nc for the normal direction on interface node c and Dsc the respective surface area. In addition, the distribution function D deals with
Table 1
Catalytic ozonation kinetic parameters for different oxidant concentrations.
ð12Þ
The species mass concentration is computed from the overall
species balance using Eq. (13) as we need to solve ns 1 transport
equations for a mixture involving ns chemical species.
C gO3 =gO3 (Nm3)
k1 (min1)
k2 (min1)
k3 (min1)
20
40
88
0.2580
0.6889
0.7413
0.2860
0.5551
0.6787
0.0059
0.0067
0.0128
0
0
0
Please cite this article in press as: R.J.G. Lopes, R.M. Quinta-Ferreira, A hybrid CFD framework for ﬂuidized bed ozonation reactors coupling interface tracking and discrete particle methods, Chem. Eng. J. (2012), http://dx.doi.org/10.1016/j.cej.2012.06.143
R.J.G. Lopes, R.M. Quinta-Ferreira / Chemical Engineering Journal xxx (2012) xxx–xxx
the numerical estimation of the Dirac-function normalized to the cell
volume. The computation of Eq. (17) was performed employing a
standard second-order ﬁnite difference method for the spatial
derivatives. Accordingly, the overall algebraic equation system is
solved once again with the application of the incomplete Cholesky
conjugate gradient algorithm. Having used a linear implicit treatment of the effective drag and added mass forces, the numerical
solution of Eq. (9) was carried out by implementing a ﬁrst-order
integration approach to track the individual motion of the catalyst
particle. The method developed by Hoomans et al. [31] has been
used to account for the non-ideal particle–particle and particle–
wall collisions.
4.2. Grid independency
The IT-DP multiphase model has been previously optimized
regarding the mesh aperture to confer grid independent results.
The numerical validation and veriﬁcation of this ﬁnite volume
method should encompass not only the optimization of mesh size
but also the accuracy obtained with different numerical schemes.
First, the scalar transport equations have been solved at different
case scenarios by evaluating the hydrodynamic variables that govern the transport phenomena occurring within the ozonation ﬂuidized bed reactor. The effect of mesh size on IT-DP simulation
results was investigated by simulating the injection of ozone bubble with 1, 2, 3, and 4 mm diameter rising in quiescent liquid in a
125 125 1000 mm three-dimensional vertical column. The IT
model replicates the injection of the ozone bubbles by generating
them immediately above the bottom plate with a speciﬁed velocity
and a period equal to the time step. Different meshes comprising
4.1 105, 8.2 105, 1.4 106, and 2.5 106 computational cells
were evaluated for obtaining the bubble diameter and bubble
velocity proﬁles, as well as the bubble detachment time. Speciﬁcally, the results attained with ﬁner meshes (1.4 106, 2.5 106)
have revealed an asymptotic trend while the coarser grids
(4.1 105, 8.2 105) underestimated the bubble velocity proﬁle.
Both IT-DP simulations have exhibited roughly identical results
for the bubble velocity, detachment time and diameter of ozone
bubbles when comparing the numerical predictions obtained with
1.4 106 and 2.5 106 cells. Conversely, the hydrodynamic data
were underpredicted when using the coarser grids with 4.1 105
and 8.2 105 cells, which pointed out that the better accuracy
can be achieved with the ﬁnest meshes without compromising
the computational expensiveness. Nevertheless the bubble behavior might be reasonably computed with 2.5 106 cells, the remainder of the VOF simulations was carried out with 1.4 106 cells
mainly due to the fact that the accuracy was not noticeably improved for such case scenario with the subsequent growth in the
number of mesh elements. Consequently, the numerically-optimized mesh with 1.4 106 has been used to investigate the inﬂuence of operating conditions and physicochemical variables on the
catalytic ozonation of phenol-like pollutants within the ﬂuidized
bed reactor.
5
Germany) that operates based on the corona effect. The ozone
gas concentration was monitored by means of an ozone gas analyzer (BMT 963 vent, BMT, Berlin, Germany). The ozone ﬂow rate
was typically 500 cm3/min with an inlet ozone concentration of
20 g/Nm3. The injection of ozone bubbles underneath the bottom
plate was accomplished in such way that none of the bubbles enters the column at the same time. This method has been used to
avoid pressure ﬂuctuations at the top of the bubble column reactor
as well as the artiﬁcial pulsing ﬂow generated by incoming bubbles, alternatively to the simultaneous introduction of ozone bubbles in all sparger holes. The surface tension has been modiﬁed
during the experiments by using different water/glycerol mixtures.
The glycerine was then diluted after each measurement with
demineralized water to obtain the ﬂuid with a viscosity approximately two times lower.
Liquid samples were withdrawn during the reaction and further
analyzed in term of total organic carbon concentration. High-performance liquid chromatography was used to measure the concentrations of the individual compounds of the model efﬂuent, as well
as some of the reaction intermediates. The samples were injected
via autosampler (Knauer Smartline Autosampler 3800), the mobile
phase (20% of methanol in water slightly acidiﬁed) was pumped
by a Knauer WellChrom K-1001 pump at a ﬂow rate of 1 mL/min
through an Eurokat H column at 85 °C, and detection was performed at 210 nm. TOC was measured with a Shimadzu 5000 Analyzer based on the combustion/nondispersive infrared gas analysis
method. Each sample was run in triplicate in order to minimize the
experimental errors. The deviation between the same sample runs
was always lower than 2% for TOC.
Ozone with different superﬁcial gas velocities in the range 0.1–
10 cm/s was fed into the reactor by means of a perforate plate. The
sparger conﬁguration has 60 holes with a diameter of 1 mm being
located in the center of the plate at square pitch of 5 mm and it is
located at the bottom of the column. The nozzles are arranged in
six equal groups of 10 needles. The sparger holes have been simulated in such a way that ozone bubbles with speciﬁc size enter the
column with a ﬁxed velocity. In order to avoid redundant binary
collisions, [0.2–2] mm bubbles were considered for the ozonation
of phenolic wastewaters and the typical distance between contiguous injections was speciﬁed to threefold of the bubble radius. A
high speed camera (Imager Pro HS CMOS camera) has been used
to record the trajectory of the ozone bubbles with an image sampling frequency at 100 s1. A 350 W halogen lamp was installed
to providing the necessary lighting conditions. The ﬁeld of view
was 12.5 21 mm (width height). The methodology to determine the bubble dimensions and its mass center involved three
positions in the ﬁeld of view: after entrance, in the middle, and before it departs. The representative bubble diameter was calculated
from the measured horizontal and vertical bubble diameter:
2
d ¼ ðdhoriz dv ert Þ1=3 .
6. Results and discussion
6.1. Effect of physicochemical properties
5. Experimental
The synthetic phenolic wastewater was prepared using
100 ppm of each of the six phenolic acids (obtained from Sigma–
Aldrich) corresponding to TOC0 = 370 mg of C/L, TPh0 (total phenolic acid content) = 350 mg of gallic acid/L, and COD0 (chemical oxygen demand) = 970 mg of O2/L. The reaction experiments have
been performed in the quasi-homogeneous ﬂow regime and the
gas oxidant was continuously fed to the reaction medium. Ozone
was generated in situ by means of pure oxygen (99.999%, Praxair,
Porto, Portugal) in an ozone generator (802 N, BMT, Berlin,
Both the surface tension and ﬂuid viscosity were investigated
on how they affect the hydrodynamic parameters of the ﬂuidized
bed reactor. In this ambit, the ozone bubble diameter and detachment time proﬁles were evaluated at different surface tension and
liquid viscosity conditions. First, the effect of surface tension on the
hydrodynamic parameters of the ozonation reactor was examined
by evaluating how the bubble diameter behaves with the increase
of superﬁcial gas velocity at different surface tensions. At
[O3] = 20 g/Nm3, T = 20 °C, and P = 1 atm, Fig. 1a shows the timeaveraged bubble diameter proﬁle at d = 4.1 102, 5.2 102,
Please cite this article in press as: R.J.G. Lopes, R.M. Quinta-Ferreira, A hybrid CFD framework for ﬂuidized bed ozonation reactors coupling interface tracking and discrete particle methods, Chem. Eng. J. (2012), http://dx.doi.org/10.1016/j.cej.2012.06.143
6
R.J.G. Lopes, R.M. Quinta-Ferreira / Chemical Engineering Journal xxx (2012) xxx–xxx
6.7 102, and 7.2 102 N/m. As the superﬁcial gas velocity increases, the ozone bubble diameter also increases at constant surface tension, and the higher the surface tension was, the larger the
ozone bubble was at constant superﬁcial gas velocity. This behavior is in agreement with the experimental evidence so different
surface tensions generate considerable distinct phenomena concerning the bubble formation, detachment time, and interstitial
ﬂow patterns. In fact, at the lowest superﬁcial gas velocity, the
bubble diameter increased as follows: 1.7, 2.2, 2.4, and 2.4 mm,
whereas the speciﬁc enlargement of ozone bubbles at the highest
superﬁcial velocity was: 4.7, 5.8, 6.5, and 6.7 mm. As expected,
the higher frequencies achieved at these conditions are conferring
the higher values of surface tension, which corresponds to the lower ozone bubble diameters.
The detachment times are depicted in Fig. 1b with different gas
superﬁcial velocities when increasing the surface tension from
d = 4.1 102 to 7.2 102 N/m. As can be seen, for the average
diameter and detachment time of ozone bubbles, both hydrodynamic variables increased gradually with the surface tension. The
detachment time delay was responsible for more gas being fed into
the growing ozone bubble as previously reported in the literature
[33]. Likewise, the increase of surface tension affected negatively
the ozone bubble generation due to long average cycle as well as
the low frequency of ozone bubble formation. Indeed, at very lower
superﬁcial gas velocities as the bubble dimensions are mainly dictated by the momentum balance including the surface tension effect, the average diameter of ozone bubbles was comparatively
high accordingly. One should bear in mind that for different viscous mixtures, the effect of surface tension is not predominant
and became negligible if one increases the gas velocity. These
mechanisms govern differently the ozone bubble size when using
low viscosity ﬂuids and the IT-DP computations for the bubble
diameter have properly mimicked those events for high surface
tension and low surface tension conditions, see Fig. 1ab. Moreover,
as to account for the buoyancy and surface tension effects, Eö was
higher than one thus emphasizing the inﬂuence of buoyancy on the
computed and experimental detachment times of ozone bubbles.
Additionally, several IT-DP simulations have been performed to
evaluate the effect of liquid viscosity on the bubble diameter and on
the bubble detachment time. At [O3] = 20 g/Nm3, T = 20 °C, and
P = 1 atm, Fig. 2a shows the time-averaged bubble diameter proﬁle
for different inlet bubble diameters. As can be seen, the larger the
inlet ozone bubble was, the faster the bubble diameter increased
at constant liquid viscosity. Indeed, at the lowest liquid viscosity
conditions (l = 0.0038 Pa s), the bubble diameter increased as follows: 4.2, 4.4, 4.5, and 5.1 mm for d0 = 0.2, 0.6, 1.0, and 1.5 mm,
a
6.2. Effect of ozone velocity
Aiming to evaluate the inﬂuence of gas velocity on the hydrodynamic operation of the ﬂuidized bed reactor, the IT-DP model has
been subsequently used to perform further simulations at different
gas superﬁcial velocities in the quasi-homogeneous and heterogeneous ﬂow regimes. Two different low gas superﬁcial velocities, uG
= 0.5 and 1 cm/s, and two high gas superﬁcial velocities, uG = 5 and
10 cm/s were evaluated on how they affect the ozone bubble and
liquid velocity ﬁeld at T = 20 °C, and P = 1 atm. Fig. 3a shows
the time-averaged radial bubble velocity proﬁle obtained at
uG = 0.5 cm/s by using the drag coefﬁcient directly computed with
IT-DP framework and applying the correlation proposed by Roghair
[27] embedded within the Euler–Lagrange model developed elsewhere [32]. In this approach, the drag coefﬁcient for single bubbles
rising in a quiescent liquid for 10 6 Re 6 1500 is expressed by
Eq. (18):
C D;1 ¼
qﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃ
€2
C D;1 Re2 þ C D;1 Eo
ð18Þ
with:
16
2
1þ
Re
1 þ 16Re1 þ 3:32Re1=2
€
4Eo
€Þ ¼
C D;1 ðEo
€ þ 9:5
Eo
C D;1 ðReÞ ¼
ð19Þ
ð20Þ
b 0.25
8
4.1e-2 N/m
5.2e-2 N/m
6.7e-2 N/m
7.2e-2 N/m
6
0.20
0.15
td, s
d B, mm
respectively. In comparison to l = 0.016 Pa s, the ozone bubble
diameters have become signiﬁcantly higher than those obtained
for low-viscosity ﬂuid conditions: 4.5, 5.0, 5.1, and 5.9 for d0 = 0.2,
0.6, 1.0, and 1.5 mm, respectively. Moreover, the highest increase
in the ozone bubble diameter was attained for the d0 = 1.5 mm, giving rise to the fact that larger bubbles are intrinsically computed
and experimentally identiﬁed for high-viscosity ﬂuids. Several ITDP simulations have been carried out to compute the detachment
time proﬁles for different viscous conditions when varying the
superﬁcial gas velocity. Fig. 2b shows the effect of liquid viscosity
on the bubble detachment time at [O3] = 20 g/Nm3, T = 20 °C, and
P = 1 atm. As can be seen, the detachment time increased considerably for constant superﬁcial gas velocity as the ﬂuid viscosity increased from l = 0.0038 to 0.016 Pa s, while the longer
detachment times were observed by decreasing the superﬁcial
gas velocity from uG = 0.5 to 0.3 m/s, and further to 0.1 m/s. From
the comparison between numerical proﬁles illustrated in Fig. 2a
and b, the time-averaged ozone bubble diameter and bubble
detachment time decreased for lower viscosity liquids.
4
0.10
2
uG=0.5m/s
0.05
uG=0.3m/s
uG=0.1m/s
0
0.0
0.2
0.4
0.6
uG, m/s
0.8
1.0
0.00
0.03
0.04
0.05
0.06
0.07
0.08
δ, N/m
Fig. 1. IT-DP predictions of (a) time-averaged bubble diameter and (b) detachment time proﬁles calculated with different surface tension ﬂuid conditions at [O3] = 20 g/Nm3,
T = 20 °C, and P = 1 atm.
Please cite this article in press as: R.J.G. Lopes, R.M. Quinta-Ferreira, A hybrid CFD framework for ﬂuidized bed ozonation reactors coupling interface tracking and discrete particle methods, Chem. Eng. J. (2012), http://dx.doi.org/10.1016/j.cej.2012.06.143
7
R.J.G. Lopes, R.M. Quinta-Ferreira / Chemical Engineering Journal xxx (2012) xxx–xxx
a
b
7
0.6
uG=0.5m/s
d0=0.2mm
0.5
d0=0.6mm
d0=1.0mm
d0=1.5mm
uG=0.1m/s
0.4
td, s
d B, mm
6
uG=0.3m/s
5
0.3
0.2
4
0.1
3
0.002
0.004
0.006
0.008
0.010
0.012
0.014
0.016
0.018
μ L, Pa.s
0.0
0.002
0.004
0.006
0.008
0.010
0.012
0.014
0.016
0.018
μ L, Pa.s
Fig. 2. IT-DP predictions of (a) Time-averaged bubble diameter and (b) detachment time proﬁles predicted by the VOF model with different liquid viscosities at [O3] = 20 g/
Nm3, T = 20 °C, and P = 1 atm.
The local gas fraction a and the drag coefﬁcient for a bubble rising in a bubble swarm was computed as follows:
CD
€1 Þð1 aÞ;
¼ ð1 þ a18Eo
C D;1
€ 6 5;
1 6 Eo
a 6 0:45
ð21Þ
As can be seen, both computed results agreed well with experimental data with a mean relative error of 11% taking into account
the low interaction regime that characterizes the multiphase ﬂow
within the ﬂuidized bed reactor. The time-averaged liquid velocity
ﬁeld is shown in Fig. 3b at the same operating conditions. The higher velocities were identiﬁed nearby the column center for both
computed results. According to these characteristic radial ﬂow patterns, one can identify the ascending motion of the liquid phase in
the center while recirculating downwards alongside the wall. The
higher gas volume fractions were attained in the center and the
lower gas volume fractions were identiﬁed near the wall. Fig. 3b
reproduced a quasi-axisymmetric pattern regarding the computed
time-averaged liquid radial velocities avoiding the artiﬁcial ﬂat and
nearly gradient-absence proﬁles for a non-optimized multiphase
reactive ﬂow model.
Fig. 4a illustrates the time-averaged radial bubble velocity proﬁles increasing the gas superﬁcial velocity to uG = 1 cm/s for the
two drag coefﬁcient frameworks at T = 20 °C, and P = 1 atm. The
computed results attained with the Roghair model resemble quite
accurately the experimental data, whereas the IT-DP model overestimated slightly the ozone bubble velocity ﬁeld with a mean relative error of 6%. Fig. 4b depicts the radial liquid velocity proﬁle
with both drag formulation approaches at the same operating conditions. Comparatively to the IT-DP model, the Roghair model has
computed the lower liquid velocities in the vertical symmetry axis
of the ﬂuidized bed notwithstanding the downward liquid velocities were almost identical near the column wall. However, both
numerical approaches have effectively reproduced the recirculation patterns with a qualitatively-featured gulf-stream proﬁle.
The IT-DP model was further used to compute bubble and liquid
velocities at high-interaction regimes by increasing the ozone
superﬁcial velocity. In Fig. 5a, the time-averaged radial bubble
velocity proﬁles were obtained at uG = 5 cm/s, T = 20 °C, and
P = 1 atm, whereas the liquid velocity proﬁles are illustrated in
Fig. 5b for the two drag coefﬁcient methodologies. The overestimation of radial bubble velocities exhibited by the IT-DP model was
almost negligible with a mean relative error of 13% both at the wall
region and in the center region of the ﬂuidized bed reactor. The
Roghair model handled agreeably the experimental data by providing conﬁdent predictions of radial bubble and liquid velocity ﬁelds,
see Fig. 5a and b. Additional CFD simulations were performed for
calculating the time-averaged radial bubble velocity proﬁles as
shown in Fig. 6a by increasing the gas superﬁcial velocity to
Fig. 3. Time-averaged radial (a) bubble and (b) liquid velocity proﬁles computed by the IT-DP CFD framework and Roghair drag coefﬁcient model at uG = 0.5 cm/s ([O3] = 20 g/
Nm3, T = 20 °C, and P = 1 atm).
Please cite this article in press as: R.J.G. Lopes, R.M. Quinta-Ferreira, A hybrid CFD framework for ﬂuidized bed ozonation reactors coupling interface tracking and discrete particle methods, Chem. Eng. J. (2012), http://dx.doi.org/10.1016/j.cej.2012.06.143
8
R.J.G. Lopes, R.M. Quinta-Ferreira / Chemical Engineering Journal xxx (2012) xxx–xxx
Fig. 4. Time-averaged radial (a) bubble and (b) liquid velocity proﬁles computed by by the IT-DP CFD framework and Roghair drag coefﬁcient model at uG = 1 cm/s
([O3] = 20 g/Nm3, T = 20 °C, and P = 1 atm).
Fig. 5. Time-averaged radial (a) bubble and (b) liquid velocity proﬁles computed by by the IT-DP CFD framework and Roghair drag coefﬁcient model at uG = 5 cm/s
([O3] = 20 g/Nm3, T = 20 °C, and P = 1 atm).
Fig. 6. Time-averaged radial (a) bubble and (b) liquid velocity proﬁles computed by by the IT-DP CFD framework and Roghair drag coefﬁcient model at uG = 10 cm/s
([O3] = 20 g/Nm3, T = 20 °C, and P = 1 atm).
uG = 10 cm/s, whereas Fig. 6b depicts the computed results for the
time-averaged radial liquid velocity proﬁles both attained at
T = 20 °C, and P = 1 atm. Consistently, the numerical results
obtained with the IT-DP model resemble accurately the experimental data for the ozone bubble velocity pattern similarly to
the application of Roghair drag coefﬁcient approach.
Please cite this article in press as: R.J.G. Lopes, R.M. Quinta-Ferreira, A hybrid CFD framework for ﬂuidized bed ozonation reactors coupling interface tracking and discrete particle methods, Chem. Eng. J. (2012), http://dx.doi.org/10.1016/j.cej.2012.06.143
9
R.J.G. Lopes, R.M. Quinta-Ferreira / Chemical Engineering Journal xxx (2012) xxx–xxx
100
10
1.6e-2 Pa.s
1.4e-2 Pa.s
5.5e-3 Pa.s
3.8e-3 Pa.s
EXP
98
C(O 3)/C(O3,0), %
TOC conversion, %
8
6
4
7.2e-2 N/m
6.7e-2 N/m
5.2e-2 N/m
4.1e-2 N/m
EXP
2
0
96
94
92
90
88
86
0
1
2
3
4
0
1
2
3
4
Time, s
Time, s
Fig. 7. Normalized total organic carbon concentration proﬁles for different surface
tension conditions at [O3] = 20 g/Nm3, T = 20 °C, and P = 1 atm.
Fig. 10. Normalized ozone concentration proﬁles for different liquid viscosity
conditions at [O3] = 20 g/Nm3, T = 20 °C, and P = 1 atm.
6.3. Detoxiﬁcation of liquid pollutants
10
TOC conversion, %
8
6
4
1.6e-2 Pa.s
1.4e-2 Pa.s
5.5e-3 Pa.s
3.8e-3 Pa.s
EXP
2
0
0
1
2
3
4
Time, s
Fig. 8. Normalized total organic carbon concentration proﬁles for different liquid
viscosity conditions at [O3] = 20 g/Nm3, T = 20 °C, and P = 1 atm.
100
C(O 3)/C(O3,0), %
98
7.2e-2 N/m
6.7e-2 N/m
5.2e-2 N/m
4.1e-2 N/m
EXP
96
94
92
90
88
86
0
1
2
3
4
Time, s
Fig. 9. Normalized ozone concentration proﬁles for different surface tension
conditions at [O3] = 20 g/Nm3, T = 20 °C, and P = 1 atm.
Under reactive ﬂow conditions, the IT-DP model has been subsequently used to carry out comprehensive simulations for evaluating the ozonation efﬁciency of liquid pollutants and how it can
be correlated with the above-presented hydrodynamic analysis.
The normalized total organic carbon concentration proﬁles are
shown in Fig. 7 for different surface tension conditions at
[O3] = 20 g/Nm3, T = 20 °C, and P = 1 atm. During the ﬁrst second
of reaction time, the computed TOC conversions were 6.4%, 6.3%,
6.0%, and 5.5% for d = 4.1 102, 5.2 102, 6.7 102, and
7.2 102 N/m. As the surface tension increases, the ozone bubble
becomes larger at constant superﬁcial gas velocity, thereby the
mass transfer area is smaller for the same ozone volumetric feed
rate. Accordingly, the mineralization rates were signiﬁcantly higher after 3 s of reaction time as follows: 7.8%, 7.5%, 6.9%, and 6.4% for
d = 4.1 102, 5.2 102, 6.7 102, and 7.2 102 N/m, respectively. The inﬂuence of liquid viscosity on the total organic carbon
removal is depicted in Fig. 8 at [O3] = 20 g/Nm3, T = 20 °C, and
P = 1 atm. As can be seen, the computed TOC conversions were
higher when low-viscosity conditions were used for the ozonation
of organic pollutants. After 1 s of reaction time, the IT-DP model
predicted the following conversions: 5.5%, 6.0%, 6.4%, and 6.9%
for l = 0.016, 0.014, 0.0055, and 0.0038 Pa s, respectively. In addition, the computed results obtained after 3 s of ozonation have revealed higher mineralization conversions: 6.5%, 7.0%, 7.7%, and
8.4%. According to the time-averaged bubble diameter proﬁles
shown in Fig. 2a and b, as the liquid viscosity increases, the ozone
bubbles become larger deteriorating the ozone mass transfer and
subsequently the reaction rates. Moreover, the highest increase
in the ozone bubble diameter was obtained for the largest inlet
diameter (d0 = 1.5 mm), thus underlining the inﬂuence exerted by
larger ozone bubbles.
Both computed and experimentally validated results regarding
the normalized ozone concentration proﬁles have been plotted in
Fig. 9 as for different surface tension at [O3] = 20 g/Nm3, T = 20 °C,
and P = 1 atm. As the surface tension increases from d = 4.1 102 to 7.2 102 N/m, the normalized ozone concentration also
increases from 91.2% to 91.7% after 1 s of reaction time, respectively. Bearing in mind the higher solubility of ozone in the liquid
phase, the mass transfer of ozone was remarkably enhanced by
conferring a distinctive time scale in comparison to pollutant concentration data depicted in Fig. 7 and 8. Concomitantly, the inﬂuence of liquid viscosity has been investigated by computing the
Please cite this article in press as: R.J.G. Lopes, R.M. Quinta-Ferreira, A hybrid CFD framework for ﬂuidized bed ozonation reactors coupling interface tracking and discrete particle methods, Chem. Eng. J. (2012), http://dx.doi.org/10.1016/j.cej.2012.06.143
10
R.J.G. Lopes, R.M. Quinta-Ferreira / Chemical Engineering Journal xxx (2012) xxx–xxx
Fig. 11. Interstitial ﬂow snapshots of ozone bubble development from IT-DP
simulations at uG = 0.1 cm/s, [O3] = 20 g/Nm3, T = 20 °C, and P = 1 atm.
normalized ozone concentration at [O3] = 20 g/Nm3, T = 20 °C, and
P = 1 atm as shown in Fig. 10. During the ﬁrst 3 s, the computed
ozone concentration comparatively decreased from 89.5% to
88.4% when the liquid viscosity decreased from l = 0.016 to
0.0038 Pa s by reinforcing the effect of this physicochemical variable on the hydrodynamic performance of the ﬂuidized bed reactor, see Fig. 2a and 2b.
6.4. IT-DP mappings of interstitial ozonation ﬂow
One should bear in mind that not only the hydrodynamic integral variables affect the overall reactor performance, but also the
mineralization parameters play a prominent effect on the feasibility of the detoxiﬁcation process. Accordingly, the inﬂuence of inlet
ozone concentration on the detoxiﬁcation efﬁciency was investigated in the range of 20–95 g/Nm3. The transient proﬁles of the total organic carbon conversion have shown that the degradation of
phenol-like compound followed a similar proﬁle as the one we obtained for different inlet gas velocities. Bearing in mind the intrinsic ozonation mechanism, the liquid pollutants have been
identiﬁed to be quickly oxidized by ozone as this is extremely reactive with compounds comprising high electronic density sites. The
catalytic ozonation started rapidly and developed a steep total organic carbon concentration proﬁle during the ﬁrst seconds of reaction time. However, as the reaction progresses the refractory
compounds were slowly oxidized revealing lower detoxiﬁcation
rates. This fact can be also ascribed to the different ozonation concentration levels attained at the ozone bubble interface as depicted
in Fig. 11. As can be seen from the development proﬁle of ozone
bubbles, the reaction system revealed to be considerably dependent on the ozone load, so the overall performance of the ﬂuidized
bed reactor can be affected for such a dynamic transport phenomena under reactive ﬂow conditions.
Additionally, the IT-DP model was used to evaluate the interstitial ﬂow patterns. In fact, the ozonation of liquid pollutants is directly related with the strength of ozone mass transfer so we aim
to analyze different computational mappings of total organic carbon concentration attained at different reaction times. The temporal evolution of total organic carbon conversion at different
operating times is shown in Fig. 12 by representing instantaneous
snapshots when the catalytic ozonation of phenolic wastewaters
was simulated at T = 20 °C and P = 1 atm. These computed results
exhibited distinctive organic matter detoxiﬁcation levels which
are inherently attributed to the ozone concentration ﬁelds depicted in Fig. 11. According to Fig. 12, the IT-DP multiphase model
also demonstrated a rapid intensiﬁcation in ozone mass transfer at
the bottom of the bubble column reactor mainly due to the higher
degree of turbulence. Indeed, the local recirculatory of liquid ﬂow
and meandering bubble plume-like structures became more intense thereby enhancing the overall decontamination efﬁciency
of phenol-like pollutants by means of catalytic ozonation within
the gas–liquid–solid ﬂuidized bed reactor.
7. Conclusions
The catalytic ozonation of liquid pollutants has been comprehensively investigated by means of an interface tracking-discrete
particle model. The hybrid CFD model was ﬁrstly developed for
Fig. 12. Interstitial ﬂow snapshots of interstitial ﬂow patterns colored by ozone concentration resulting from IT-DP simulations for bubble plume structures at uG = 0.1 cm/s,
[O3] = 20 g/Nm3, T = 20 °C, and P = 1 atm).
Please cite this article in press as: R.J.G. Lopes, R.M. Quinta-Ferreira, A hybrid CFD framework for ﬂuidized bed ozonation reactors coupling interface tracking and discrete particle methods, Chem. Eng. J. (2012), http://dx.doi.org/10.1016/j.cej.2012.06.143
R.J.G. Lopes, R.M. Quinta-Ferreira / Chemical Engineering Journal xxx (2012) xxx–xxx
the gas–liquid–solid reaction system to query the inﬂuence of
physicochemical variables on the hydrodynamic operation by
addressing different surface tension and liquid viscosity conditions. Afterwards, several numerical simulations have been accomplished under distinct ﬂow regimes for evaluating the effect of
ozone velocity on the time-averaged bubble diameter and detachment time. Here, the increase of ﬂuid viscosity has been found to
generate larger ozone bubbles thus deteriorating the ozone mass
transfer, and subsequently the mineralization rates. Under different process conditions, the IT-DP handled agreeably the experimental data both for hydrodynamic and reactive ﬂow
parameters. Finally, the interstitial ﬂow mappings of the gas–liquid–solid reactor gave rise to prominent ozone concentration
ﬁelds by unveiling a heterogeneous degree of transport phenomena. This fact has been ascribed to the bubble plume hydrodynamic
structures thereby affecting the mineralization efﬁciency attained
by the catalytic ozonation reactor.
Acknowledgment
The authors gratefully acknowledged the ﬁnancial support of
Fundação para a Ciência e Tecnologia, Portugal.
References
[1] C.N. Satterﬁeld, G.A. Huff, Effects of mass transfer on Fischer–Tropsch synthesis
in slurry reactors, Chemical Engineering Science 35 (1980) 195–202.
[2] A.E. Cover, W.C. Schreiner, G.T. Skaperdas, Kellogg Coal Gasiﬁcation Process
Chemical Engineering Progress 69 (1973) 31–36.
[3] M.P. Dudukovic, F. Larachi, P.L. Mills, Multiphase reactors – revisited, Chemical
Engineering Science 54 (1999) 1975–1995.
[4] L.-S. Fan, Gas–Liquid–Solid Fluidization Engineering, Butterworth Series in
Chemical Engineering, Butterworth, Stoneham, MA, 1989.
[5] E. Delnoij, J.A.M. Kuipers, W.P.M. van Swaaij, Computational ﬂuid dynamics
applied to gas–liquid contactors, Chemical Engineering Science 52 (1997)
3623–3638.
[6] R. Krishna, J.M. van Baten, Scaling up bubble column reactors with the aid of
CFD, Chemical Engineering Research and Design 79 (2001) 283–309.
[7] A.A. Kulkarni, J.B. Joshi, Bubble formation and bubble rise velocity in gas–liquid
systems: a review, Industrial Engineering Chemistry Research 44 (2005) 5873–
5931.
[8] Y. Pan, M.P. Dudukovic, Numerical investigation of gas-driven ﬂow in 2-D
bubble columns, AICh.E. Journal 46 (2000) 434–449.
[9] D. Pﬂeger, S. Becker, Modeling and simulation of the dynamic ﬂow behavior in
a bubble column, Chemical Engineering Science 56 (2001) 1737–1747.
[10] E. Delnoij, J.A.M. Kuipers, W.P.M. Van Swaaij, A three-dimensional CFD model
for gas–liquid bubble columns, Chemical Engineering Science 54 (1999) 2217–
2226.
[11] S. Lain, D. Bröder, M. Sommerfeld, Experimental and numerical studies of the
hydrodynamics in a bubble column, Chemical Engineering Science 54 (1999)
4913–4920.
[12] R.O. Fox, On the relationship between Lagrangian micromixing models and
computational ﬂuid dynamics, Chemical Engineering Progress 37 (1998) 521–
535.
11
[13] W. Gerlinger, K. Schneider, L. Falk, H. Bockhorn, Numerical simulation of the
mixing of passive and reactive scalars in two-dimensional ﬂows dominated by
coherent vortices, Chemical Engineering Science 55 (2000) 4255–4269.
[14] G.Q. Yang, B. Du, L.S. Fan, Bubble formation and dynamics in gas–liquid–solid
ﬂuidization – a review, Chemical Engineering Science 62 (2007) 2–27.
[15] J.U. Brackbill, D.B. Kothe, C.A. Zemach, A continuum method for modeling
surface tension, Journal of Computational Physics 100 (1992) 335–354.
[16] C.W. Hirt, B.D. Nichols, Volume of ﬂuid (VOF) method for the dynamics of free
boundaries, Journal of Computational Physics 39 (1981) 201–225.
[17] S. Osher, J.A. Sethian, Front propagating with curvature-dependent speed:
algorithms based on Hamilton–Jacobi formulations, Journal of Computational
Physics 79 (1988) 12–49.
[18] Y. Li, J. Zhang, L.-S. Fan, Numerical simulation of gas–liquid–solid ﬂuidization
systems using a combined CFD–VOF–DPM method: bubble wake behaviour,
Chemical Engineering Science 54 (1999) 5101–5107.
[19] N. Yang, W. Wang, W. Ge, J. Li, Modeling of meso-scale structures in
particle–ﬂuid systems: the EMMS/CFD approach, China Particuology 3
(2005) 78–79.
[20] A. Renzo, F. Maio, Homogeneous and bubbling ﬂuidization regimes in DEM–
CFD simulations: Hydrodynamic stability of gas and liquid ﬂuidized beds,
Chemical Engineering Science 62 (2007) 116–130.
[21] L. Mazzei, P. Lettieri, CFD simulations of expanding/contracting homogeneous
ﬂuidized beds and their transition to bubbling, Chemical Engineering Science
63 (2008) 5831–5847.
[22] R. Panneerselvam, S. Savithri, G.D. Surender, CFD simulation of hydrodynamics
of gas–liquid–solid ﬂuidised bed reactor, Chemical Engineering Science 64
(2009) 1119–1135.
[23] L. Montastruc, J.P. Brienne, I. Nikov, Modeling of residence time distribution:
application to a three-phase inverse ﬂuidized bed based on a Mellin transform,
Chemical Engineering Journal 148 (2009) 139–144.
[24] K. Sivaguru, K.M. Meera Sheriffa Begum, N. Anantharaman, Hydrodynamic
studies on three-phase ﬂuidized bed using CFD analysis, Chemical Engineering
Journal 155 (2009) 207–214.
[25] Y. Behjat, S. Shahhosseini, M. Marvast, Investigation of catalyst particle
hydrodynamic and heat transfer in three phase ﬂow circulating ﬂuidized
bed, International Communications in Heat and Mass Transfer 38 (2011)
100–109.
[26] X. Huang, CFD modeling of liquid–solid ﬂuidization: effect of drag correlation
and added mass force, Particuology 9 (2011) 441–445.
[27] I. Roghair, Y.M. Lau, N.G. Deen, H.M. Slagter, M.W. Baltussen, M. Van Sint
Annaland, J.A.M. Kuipers, On the drag force of bubbles in bubble swarms at
intermediate and high Reynolds numbers, Chemical Engineering Science 66
(2011) 3204–3211.
[28] S.O. Unverdi, G. Tryggvason, A front-tracking method for viscous,
incompressible multi-ﬂuid ﬂows, Journal of Computer Physics 100 (1992)
25–37.
[29] A. Prosperetti, Navier–Stokes numerical algorithms for free-surface ﬂow
computations: an overview, Drop Surface Interactions 21 (2001).
[30] E. Delnoij, F.A. Lammers, J.A.M. Kuipers, W.P.M. van Swaaij, Dynamic
simulation of dispersed gas–liquid two-phase ﬂow using a discrete bubble
model, Chemical Engineering Science 52 (9) (1997) 1429–1458.
[31] B.P.B. Hoomans, J.A.M. Kuipers, W.J. Briels, W.P.M. van Swaaij, Discrete particle
simulation of bubble and slug formation in a two-dimensional gas-ﬂuidised
bed: a hard sphere approach, Chemical Engineering Science 51 (1996) 99–108.
[32] R.J.G. Lopes, R.M. Quinta-Ferreira, Detoxiﬁcation of high-strength liquid
pollutants in an ozone bubble column reactor: gas–liquid ﬂow patterns,
interphase mass transfer and chemical depuration, Chemical Engineering
Journal 172 (2011) 476–486.
[33] D. Gerlach, G. Tomar, G. Biswas, F. Durst, Comparison of volume-of-ﬂuid
methods for surface tension-dominant two-phase ﬂows, International Journal
of Heat Mass Transfer 49 (2006) 740–754.
Please cite this article in press as: R.J.G. Lopes, R.M. Quinta-Ferreira, A hybrid CFD framework for ﬂuidized bed ozonation reactors coupling interface tracking and discrete particle methods, Chem. Eng. J. (2012), http://dx.doi.org/10.1016/j.cej.2012.06.143