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

© Copyright 2019