0263–8762/06/$30.00+0.00 # 2006 Institution of Chemical Engineers Trans IChemE, Part A, August 2006 Chemical Engineering Research and Design, 84(A8): 677– 690

www.icheme.org/cherd doi: 10.1205/cherd.05178

CFD SIMULATION OF FLOW AND AXIAL DISPERSION IN EXTERNAL LOOP AIRLIFT REACTOR S. ROY, M. T. DHOTRE and J. B. JOSHI
Institute of Chemical Technology, University of Mumbai, Matunga, Mumbai, India

I

n the present work, a CFD simulation, using two-fluid model, has been performed to obtain the flow pattern and the gas hold-up in external loop air-lift reactors (EL-ALR). An agreement was observed between the predicted and the experimental data available in the published literature. The validated CFD model has been extended for the simulation of mixing and the axial dispersion coefficient. CFD predictions of mixing time are well in agreement with the experimental values published in literature, and the results of axial dispersion coefficient also show good agreement with the experimental values published in literature except in the case of Vial et al. (2005), where the CFD model underpredicts. A systematic numerical study was undertaken to understand the relative contribution of convection and eddy diffusion to the mixing and axial dispersion. Keywords: CFD; mixing; axial dispersion coefficient; airlift reactor.

INTRODUCTION

namely the downcomer, the riser and the gas disengagement section (Figure 1). Sparging of the gas phase can be accomplished in the riser and/or in the downcomer. Sparging in the riser leads to a gas – liquid up-flow in the riser and liquid downflow in the downcomer, leading to a well-directed circulation pattern. On the other hand, sparging in downcomer results in substantial increases in gas –liquid interfacial area, absorption rate and a reduction of power consumption per unit volume of the reactor. In spite of frequent use of EL-ALR, present design procedure is highly empirical because of the complex fluid dynamics prevailing in the reactor. Coexistence of both single phase and multiphase make EL-ALR hydrodynamically most challenging, and thus modelling it as a single reactor poses a challenging problem (Joshi et al., 1990). Each section, namely riser, downcomer and disengagement section, is to be understood and modelled separately. For reliable design of EL-ALR, it is desirable to have a detailed knowledge of the flow pattern. Further, this information can be extended for estimation of various design objectives like mixing time, axial dispersion coefficient, heat transfer coefficient, and so on. Further, the axial dispersion in the riser of the EL-ALR can be varied from that of a bubble column in one extreme to a pipe flow in other extreme, which makes EL-ALR very interesting to study. To quantify this wide range of behaviour we have to first understand the hydrodynamics (profiles of mean velocity, holdup, eddy diffusivity, and so on) of the reactor. During the past 30 years, there have been continuous and vigorous attempts in the direction of reducing empiricism. In particular, developments in computational fluid dynamics (CFD) paved the way for this in the past two decades because

Over the past three decades many researchers across the world have put their quality time on understanding the hydrodynamic behaviour of airlift reactors, the most popular modifications of bubble column. The modifications were proposed to overcome two major limitations: (1) complete back-mixing in the liquid phase and (2) high pressure drop. Moreover, this class of reactor is very attractive for use in the chemical process industry and biotechnology due to their design flexibility, low power requirement, and with further advantages of good mass and heat transfer. The more remarkable advantage of these reactors is the control over liquid circulation. The well directive liquid circulation pattern obtained in these type of reactors, reduce the shear stress to a considerable amount, thus facilitating the cultivation of shear sensitive organisms, as a result these reactors are widely used in the biochemical industry and for waste water treatment. The large reactor volumes needed for the economic production of chemicals and biological organism, or for effluent treatment, can be achieved with high aspect ratio (height-to-diameter ratio) in external loop air lift reactors (EL-ALR) without sacrificing mixing (contrary to conventional bubble columns). The EL-ALR have greater flexibility (Weiland and Onken, 1981; Siegel et al., 1986; Verlaan et al., 1986; Joshi et al., 1990) and its performance can be manipulated better by controlling parameters for the individual sections
Correspondence to: Professor J. B. Joshi, Institute of Chemical Technology, University of Mumbai, Matunga, Mumbai 400 019, India. E-mail: [email protected]

677

678

ROY et al.

Figure 1. Schematic of the external loop airlift reactor. (1,5—gas inlet, 2— downcomer, 3—riser, 4—bottom bend, 6—gas disengagement section.

of the spectacular progress in the digital computing. The CFD has now emerged as a powerful tool for solving the governing equations of change over wide range of Reynolds number covering molecular, inertial and turbulent modes of transport. A balanced combination of CFD and experiments gives fairly good knowledge of flow pattern in the equipment under consideration. Several recent publications have established the potential of CFD for describing the hydrodynamics of bubble column (Jakobsen et al., 1997; Joshi, 2001; Ekambara and Joshi, 2003). In view of this, an attempt has been made in the present work to use the CFD tool for EL-ALR and validate it with the available experimental data in the literature. Further, it was thought desirable to extend the flow information for the prediction of mixing time (umix) and axial dispersion coefficient (DL) and to validate the model by comparing with the experimental data available in the published literature. The main objective of this paper is to bring out the following points: (1) relative contribution of convection and eddy diffusion to the mixing time and axial dispersion, (2) understanding the basis of an order of magnitude lower value of dispersion coefficient in EL-ALR as compared with bubble column, (3) the effect of area ratio, one of the most important design parameter on mixing time and axial dispersion.

PREVIOUS WORK A wide range of experimental studies of airlift hydrodynamics has been conducted in recent years. A large number of investigators (Hsu and Dudukovic, 1980; Merchuck and Stein, 1981; Nichol, 1984; Merchuck, 1986; Siegel et al., 1986; Verlaan et al., 1986; Popovic and Robinson, 1987a, b, 1989; Akita et al., 1988; Chisti and Moo Young, 1988; Young et al., 1991; Okada et al., 1993; Bentifraouine et al., 1997; Bendjaballah et al., 1999; Vial et al., 2005)

have published papers reporting the variation of gas holdup, liquid circulation velocity and mass transfer coefficient, primarily with VG and AD/AR for a wide range of geometry, gas velocities, gas liquid systems (viz. CMC solution, NaCl solution). Thus, the hydrodynamic quantities that have been examined to date include phase velocities, flow regimes and gas-phase volume fraction. As regards to mixing and axial dispersion, a very little information is available in the published literature on the experimental investigations in EL-ALR. Each of these quantities is influenced by the various design and operating features (like range of VG, geometrical configurations, nature of liquid phase). In the published literature confusions still persists (untill 2005) regarding the dependencies of mixing time and axial dispersion on these design and operating parameters. A brief account of these studies has been given below: Weiland and Onken (1981) studied hydrodynamics, mass transfer of oxygen from gas to liquid, axial dispersion (in the riser) and mixing behaviour in tubular loop fermenter with a overall height of over 10 m. They obtained the mixing time for 10% inhomogeneity, and found it to be a function of superficial gas velocity (VG) only. It decreases rapidly at low VG range (up to 0.02 m s21) and tend to level off with increasing VG. The Bodenstein number (Bo) was calculated by fitting the model transfer function [ratio of the laplace transformed output to the input (ideal Dirac pulse) concentration] to the experimental transfer function obtained by taking laplace transform of experimental concentration of tracer at two locations in the riser. They showed that axial dispersion coefficient (DL) decreases with increasing liquid velocity (in the range 0.07 – 0.28 m s21) for a constant VG and DL increases with increase in VG. They have also proposed a correlation for axial dispersion considering the influence of VG on axial dispersion to be comparable with that in bubble column. Bello et al. (1984) obtained the liquid circulation velocity, liquid phase mixing time and Bo for three EL-ALRs, having different AD/AR ratio (0.11, 0.25 and 0.44). Mixing time was calculated for 1% inhomogeneity. They concluded that the ratio of mixing time to circulation time is independent of VG and dependent only on AD/AR ratio. The Bo has been calculated by fitting the experimental response to the theoretical response for an open – open vessel as given in Levenspiel (1995). They reported the Bo for whole reactor and its variation with VG (in the range 0.014 – 0.086 m s21) and AD/AR ratio, which did not give any idea of the value of DL in the riser and its variation with different parameters. Verlaan et al. (1989) investigated the axial dispersion in different sections of an EL-ALR and mixing time. They have calculated mixing time for 5% inhomogeneity and observed similar kind of behaviour as that reported by Weiland and Onken (1981). For each section, the predicted outlet concentration was determined by calculating Fourier transformed transfer function and then taking inverse Fourier transform of the same. The Bo was then calculated by fitting experimental output concentration of the section with the respective predicted outlet concentration. They found DL to be almost independent of the VG (in the range 0.014– 0.14 m s21), and concluded that the turbulence induced by the liquid velocity forms the main contribution to dispersion. Kawase et al. (1994) measured mixing time, circulation time and downcomer linear liquid velocity for two EL-ALR, having different AD/AR ratio (0.204 and 0.458).

Trans IChemE, Part A, Chemical Engineering Research and Design, 2006, 84(A8): 677– 690

CFD SIMULATION OF FLOW AND AXIAL DISPERSION Mixing time has been calculated for 5% inhomogeneity. They observed that mixing time to decrease with increasing VG (in the range 0.01– 0.045 m s21). However, beyond 0.02 m s21 the mixing time was found to attain a constant value. They also reported the increase of mixing time with increase in AD/AR ratio. Gavrilescu and Tudose (1997) have performed extensive experimental studies on mixing and axial dispersion behaviour of the liquid phase in two pilot scale EL-ALR with different values of AD/AR ratio (0.040 and 0.1225). Mixing time was calculated for 5% inhomogeneity. They found that mixing time decreases with increasing VG (in the range 0.01 –0.12 m s21), and it increases with increases in AD/AR ratio. They concluded that the mixing time as a measure of macroscale mixing promoted by convective mechanism. DL for different sections of EL-ALR was calculated following Verlaan’s method. They have found an increase in DL with increase in VG (in the range 0.01– 0.12 m s21). They also reported a decreases in DL with an increase in AD/AR ratio. Their results showed an overriding contribution of convective transport on DL. Vial et al. (2005) have performed experimental and theoretical studies of axial dispersion in the riser of two different EL-ALR (AD/AR ¼ 0.25 and 0.28). Concentration profile at the outlet of the riser has been estimated by solving the transport equation for a passive scaler. The value of DL was then obtained by fitting experimental data to the estimated concentration profile. They have reported very high (two to five times) value of experimental axial dispersion coefficient in compassion to all those published to date. Further, they have concluded that the axial dispersion is mainly controlled by the liquid phase turbulence and the convection has practically no influence. A summary of the published work has been given in Table 1. From Table 1 and foregoing discussion, it is

679

clear that further investigation are required for understanding the following points: (1) wide variation of results on axial dispersion coefficients (DL), (2) the effect of area ratio and the riser diameter on the mixing time (umix) and DL, (3) the relative contribution of convection and turbulence (eddy diffusivity) on umix and DL. For achieving this objective it was thought desirable to simulate the flow pattern over wide range of EL-ALR geometry. Further, extensive comparison has been made with the published results of Weiland and Onken (1981), Bello et al. (1984), Verlaan et al. (1989), Kawase et al. (1994) Gavrilescu and Tudose (1997), Vial et al. (2005). MATHEMATICAL MODEL The equations of continuity, motion, turbulent kinetic energy (k) and turbulent energy dissipation (1) for threedimensional cylindrical co-ordinate system can be represented by the following governing equation: @ (r [k ) þ r  (rk [k uk ) ¼ 0 @t k @ (r [k uk ) þ r  (rk [k uk uk ) @t k ¼ r  ( [k tk )  [k rP þ [k rk g þ MI,k

(2)

The terms on the right hand side of equation (2) are respectively representing the stress, the pressure gradient, gravity and the ensemble averaged momentum exchange between the phases, due to interface forces. The pressure is defined equal in both phases. The stress term of phase k is described as follows:
 2 T tk ¼ meff,k ruk þ (ruk )  I(r  uk ) (3) 3

Table 1. Axial dispersion coefficient: experimental details of previous work. Authors Weiland and Onken (1981) Bello et al. (1984) Verlaan et al. (1989) Kawase et al. (1994)

Gavrilescu and Tudose (1997) Vial et al. (2005)

Range of VG

Tracer used

Measuremet technique

Air–water þ NaCl (0.1 and 1 M) soln./propanol soln./saccharose soln. (m ¼ 0.004–0.016 Pa s21) Air–water Air–water þ NaCl soln.

0.005–0.05 m s21

Aqueous NaCl solution

Conductivity measurement

0.0137–0.086 m s21

Aqueous NaCl and methylene blue dye

Air–water

0.0052–0.171 m s21

Acid þ Base

Conductivity measurement and Colorimetric measurement pH Measurement

Air–water Air–glycerine/CMC (0.05 w%)/CMC (0.15 w%)/Xanthan gum mCMC ¼ 0.0069–0.04 Pa s21 mXanthan ¼ 0.19 Pa s21 Air–water

0.005–0.045 m s21

HCl

pH measurement

ELALR-P 0.01–0.12 m s21

ELALR-P Aqueous NaCl

ELALR-P Conductivity measurement

R1 0.009– 0.25 m s21 R2 0.011– 0.113 m s21

Saturated NaCl

Conductivity measurement

Details of geometry

Details of gas–liquid system

External loop DR ¼ 0.1 m DD ¼ 0.05 m Height ¼ 10.0 m External loop DR ¼ 0.152 m DD ¼ 0.051, 0.076, 0.102 m Hdisp ¼ 1.8 External loop DR ¼ 0.2 m DD ¼ 0.1 m Height ¼ 3.23 m External loop DR ¼ 0.155 m DD ¼ 0.105, 0.07 m HR ¼ 1.37 m HD ¼ 0.60, 0.80 m Vol ¼ 0.026, 023 m3 External loop (ELALR-P)-DR ¼ 0.2 m DD ¼ 0.07, 0.04 m HR ¼ 4.70 m HD ¼ 4.40 m External loop R1 DR ¼ 0.1 m DD ¼ 0.05 m H ¼ 2.75 m External loop R2 DR ¼ 0.15 m, DD ¼ 0.08 m H ¼ 6.0 m

Air–water

Trans IChemE, Part A, Chemical Engineering Research and Design, 2006, 84(A8): 677 –690

680

ROY et al.

where, meff,k is the effective viscosity. The effective viscosity of the liquid phase is composed of three contributions: the molecular viscosity, the turbulent viscosity and an extra term due to bubble induced turbulence.

meff,L ¼ mL þ mT,L þ mBI, L

(4)

The calculation of the effective gas viscosity is based on the effective liquid viscosity (Jakobsen et al., (1997) as follows:

meff,G ¼

rG m rL eff,L

(5)

There are several models to take account of the turbulence induced by the movement of the bubbles. In this study the model proposed by Sato and Sekoguchi (1975) was used.

mBI,L ¼ rL Cm,BI [G dB juG  uL j

(6)

Thakre and Joshi (1999) have also pointed out the negligible effect of virtual mass force. The drag force is given as follows: 3 CD juG  uL j(uG  uL ) MD,L ¼  [G rL 4 dB

(12)

where CD is the drag coefficient taking into account the character of the flow around the bubble. We have used following drag coefficient (Zhang and Vanderheyden, 2002) CD ¼ 0:44 þ

24 6 pffiffiffiffiffiffi þ Re 1 þ Re

(13)

The bubble size was chosen in such a way that, it satisfies the average slip of the gas liquid system under consideration. In general dB value varies from 2 to 6 mm. The effect if shearing motion in the liquid phase on the movement of the gas phase is modelled through the lift force. MI,L ¼ CL [G rL (uG  uL )  r  uL

(14)

The turbulent eddy diffusivity is formulated as follows: 2

meff,L ¼ rL Cm

k 1

(7)

The turbulent kinetic energy k and its energy dissipation rate 1 are calculated from their governing equations: @ (r [L k) þ r  (rL [L uL k) @t L
 meff,L ¼ r  [L rk þ [L (G  rL 1) (8) sk @ (r [L 1) þ r  (rL [L uL 1) @t L
 mL;eff 1 r1 þ [L (C 11 G  C 12 rL 1) ¼ r  [L k s1 (9) The model constants are Cm ¼ 0.09; Cm,BI ¼ 0.6; sk ¼ 1.00; s1 ¼ 1.00; C11 ¼ 1.44, C12 ¼ 1.92. The term G in above equation is the production of turbulent kinetic energy and described by G ¼ t L :ruL

(10)

The total interfacial force acting between two phases may arise from several independent physical effects: MI,L ¼ MI,G ¼ MD,L þ ML,L þ MVM,L þ MTD,L

(11)

The forces indicated above respectively represent the interphase drag force, lift force, virtual mass force and turbulent dispersion force. In the present hydrodynamic model all forces (the drag force, lift force and turbulence dispersion forces) except the virtual mass force has been used. According to Hunt et al. (1987) the contribution of virtual mass force becomes negligible for diameter greater that 0.15 m. Deen et al. (2001), Sokolichin et al. (2004),

where CL is the lift coefficient which is set to 0.1. The turbulent dispersion force, derived by Lopez de Bertodano (1991), is based on the analogy with molecular movement. It approximates a turbulent diffusion of the bubbles by the liquid eddies. It is formulated as MTD,L ¼ MTD,G ¼ CTD rL kr [L

(15)

where k is the liquid turbulent kinetic energy per unit of mass. The turbulent dispersion coefficient CTD ¼ 0.1 have been used in the simulation.

Model Formulation for the Dispersion The mixing process occurs due to the transport at three levels: molecular, eddy and bulk (convection). Usually, the bulk motion (or bulk diffusion) is superimposed on either molecular or eddy diffusion or both. The convective transport can be characterized by the mean velocity components in the liquid phase. The dispersive transport can be characterized by eddy diffusivity in the liquid phase. In view of the controlling role of the convection and eddy diffusion, the solution of convection –diffusion equation for the liquid phase can give an insight into the mechanism of axial mixing and hence a rational expression for the effective axial dispersion coefficient. The transient mass balance for a tracer substance in a three-dimensional axisymmetric cylindrical co-ordinate system is given by the following equation: @ ( [L c) þ r  ( [L uL c) ¼ r  (Deff [L r  c) @t

(16)

The dispersion of the tracer is expressed by the transient mass balance [equation (16)]. The effective diffusion coefficient, Deff, is given by Deff ¼ Dm þ Dt

(17)

where Dm is the molecular diffusion coefficient and Dt is the eddy diffusion coefficient. The value of Dm was taken as 0.00005 m2/s21 in all the simulations.

Trans IChemE, Part A, Chemical Engineering Research and Design, 2006, 84(A8): 677– 690

CFD SIMULATION OF FLOW AND AXIAL DISPERSION

681

NUMERICAL DETAILS

Material Balance

The boundary conditions were defined as follows: at the gas inlet the spargers were incorporated by creating source points and specifying the mass flow rate and inlet hole velocity which was calculated by

The drift-flux model of Zuber and Findlay (1969) is given by the following equation:

Vo ¼

VG =[G ¼ C0 (VG þ VL ) þ C1

(21)

where

VG D 2 Ndo2

(18)

Along the walls, no-slip boundary conditions were adopted. At the outlet of column, outlet boundary conditions were used and atmospheric pressure was specified. The steady state equations were solved using CFX 5.7.1. The solution to the diffusion-convection equation (16) must satisfy the following initial and boundary conditions: c(r, z, u, 0) ¼ C00 @c(0, z, t) @c(R, z, t) ¼ ¼0 @r @r

(19) for all z and t

The continuity and momentum balance equations were coupled using the SIMPLEC algorithm. Three-dimensional steady state simulations were performed. Pressure was never under-relaxed, as required by the SIMPLEC algorithm. The resulting linear equation sets were solved by using suitable single-block or multi-block versions of the strongly implicit procedure for the momentum and scalartransport equations. The hybrid-upwind discretization scheme was used for the convective terms. Simulations were carried out with a tetrahedral mesh (10 000 – 50 000 nodes depending on the geometry).

C0 ¼

k [G ( [G uG þ V L )l k [G lkV G þ V L l

(22)

C1 ¼

k [G [L V S l k [G l

(23)

and

The parameters C0 and C1 are the drift– flux constants. C0 represents the hold-up profile and C1 the bubble rise velocity. The most fortunate characteristic feature of this model is that the values of C0 and C1 are practically independent of the column diameter (of course when D . 150 mm and the sparger region is exceeded). Therefore, ˆ for a given gas – liquid system, a few measurements  G with respect to VG and VL (over the range of interest) of [ in a small diameter column (150 mm) enable the estimation of C0 and C1. It is important that equation (21) holds for an extreme case of homogeneous regime having C0 ¼ 1. The set of steady state governing equations were solved numerically. Once the flow pattern in the column is obtained by the iterative procedure, the results of velocity profiles and eddy diffusivity were further used for the

Energy Balance All the predicted flow patterns must satisfy the energy balance. The rate of energy supply from the gas phase to the liquid phase is given by the following equation (Dhotre and Joshi, 2004):

p  G Vs  E1 ¼ D2 (rL  rG )gHdisp [L ½VG þ(CB 1)[ 4

(20)

When bubbles rise, the pressure energy is converted into turbulent kinetic energy. A fraction of CB is considered to get transferred to the liquid phase; the rate of energy given by equation (20) is finally dissipated in the turbulent liquid motion. It must be emphasized at this stage that, whatever may be the value of CB, the energy balance must be satisfied. When k– 1 model is used for the prediction of flow pattern, we get values of 1 (energy dissipation rate per unit mass) as one of the answers. From this 1 field, the total energy dissipation rate can be calculated by suitable volume integration. The total energy dissipation rate must equal the energy-input rate given by equation (20). The pertinent detailed discussion has been provided by Dhotre and Joshi (2004).

Figure 2. Comparison between predicted and experimental liquid circulation velocity: S, Weiland and Onken (1981); þ, Bello et al. (1984) [AD/ j , Bello et al. (1984) AR ¼ 0.11]; O, Bello et al. (1984) [AD/AR ¼ 0.25];  [AD/AR ¼ 0.44]; A, Verlaan et al. (1989); D, Gravrilescu and Tudose (1997) [AD/AR ¼ 0.12]; †, Gravrilescu and Tudose (1997) [AD/ AR ¼ 0.04]; B Vial et al. (2005) [AD/AR ¼ 0.25]; , Vial et al. (2005) [AD/AR ¼ 0.28].

Trans IChemE, Part A, Chemical Engineering Research and Design, 2006, 84(A8): 677 –690

682

ROY et al.

estimation of the axial dispersion coefficient by solving the transient mass balance equation for an inert tracer.

The liquid circulation velocity has been identified as a key design parameter for air lift reactors by many researchers (Weiland and Onken, 1981; Verlaan et al., 1989; Kochbeck et al., 1992). The performance of an EL-ALR and the reliability of its operation largely depend on the extent of liquid circulation, which is the direct manifestation of the AD/AR ratio. It governs the design parameters such as gas –liquid effective interfacial area, rate of mass transfer, extent of liquid phase mixing and the rate of heat transfer at the wall. Thus, it is very important to predict accurately the liquid circulation rate in order to have a good understanding of the operational performance and to achieve optimum design of these contactors for various applications (Bello et al., 1984). In view of this, simulations have been carried out for experimental data of Weiland and

Onken (1981), Bello et al. (1984), Verlaan et al. (1989), Gavrilescu and Tudose (1997) and Vial et al. (2005). The average liquid velocity has been predicted and comparison has been shown in the Figure 2. It can be seen that excellent parity has been obtained for all the cases. For validating the model for the flow pattern, simulations have been carried out for the experimental conditions of Young et al. (1991) and the comparison has been shown in Figures 3 and 4 and details have been given in Table 2 for two downcomer diameters of 140 mm and 89 mm. A good agreement was found between the model predictions and the experimental data for the profiles of liquid velocity as well as the gas hold-up. Some deviation, however, can be seen (Figures 3 and 4) in the near wall region. This is because the wall region has not been resolved in the simulation. The effect of downcomer diameter (in turn AD/AR) is clearly visible in the results. As the downcomer diameter increases (AD/AR increases), the axial liquid velocity also increases, owing to the less pressure drop in the downcomer. In view of this favourable comparison, it was thought desirable to extend the model for the prediction of design objectives like the mixing time and the axial dispersion coefficient.

Figure 3. Comparison of CFD predictions and experimental data of Young et al. (1991) for (a) hold-up profile, (b) axial liquid velocity for DD ¼ 8.9 cm; O, VG ¼ 0.0096 m s21; †, VG ¼ 0.021 m s21; V, VG ¼ 0.047 m s21; B, VG ¼ 0.084 m s21.

Figure 4. Comparison of CFD predictions and experimental data of Young et al. (1991) for (a) hold-up profile, (b) axial liquid velocity for DD ¼ 14 cm; O, VG ¼ 0.0096 m s21; †, VG ¼ 0.021 m s21; V, VG ¼ 0.047 m s21; B, VG ¼ 0.084 m s21.

RESULTS AND DISCUSSION Comparison of the Flow Pattern with the Experimental Data

Trans IChemE, Part A, Chemical Engineering Research and Design, 2006, 84(A8): 677– 690

CFD SIMULATION OF FLOW AND AXIAL DISPERSION

683

Table 2. Comparison of CFD predictions with the experimental data of Young et al. (1991). Material balance
VG (m s21)

AD/AR

VL (m s21)

[G

Energy balance



LHS

RHS

LHS

RHS

0.49

0.01 0.023 0.050 0.087

(0.0096) (0.021) (0.047) (0.082)

0.265 0.338 0.44 0.619

(0.266) (0.343) (0.448) (0.617)

0.024 0.049 0.080 0.121

(0.024) (0.049) (0.080) (0.120)

0.40 0.428 0.587 0.683

0.319 0.411 0.554 0.610

0.0826 0.1967 0.4067 0.7639

0.0828 0.2065 0.3591 0.7467

0.02

0.0095 0.021 0.046 0.085

(0.0096) (0.021) (0.047) (0.084)

0.136 0.20 0.286 0.330

(0.136) (0.199) (0.286) (0.344)

0.033 0.059 0.102 0.149

(0.033) (0.059) (0.102) (0.147)

0.291 0.355 0.460 0.570

0.285 0.342 0.428 0.517

0.066 0.142 0.297 0.518

0.069 0.140 0.303 0.517

The bracketed numbers indicate the experimental values.
¯ 1; LHS ¼ V¯G/[G; RHS ¼ C0(V¯G þ V¯L) þ C

LHS ¼ volume integral of energy dissipation, RHS ¼ energy input.

Comparison Between the Predicted and the Experimental Values of Mixing Time The published experimental work has been summarized in Table 1, which gives the details pertaining to column diameter, column height, gas – liquid system, range of VG, tracer and the measurement techniques employed for the measurement of axial dispersion coefficient. The

simulations have been carried out for the experimental conditions of Weiland and Onken (1981), Bello et al. (1984), Verlaan et al. (1989), Kawase et al. (1994) and Gavrilescu and Tudose (1997). The mixing time was determined directly from the tracer response curves for certain degree of inhomogeneity as reported by the different authors. Figure 5 shows a typical tracer response curve at the outlet of riser. The comparison between the predicted results and experimental data is given in Tables 3– 7. The parity of the mixing time has been shown in Figure 6. It can be seen that the predictions are in good agreement with the experimental data. At this stage we would like to bring out one important point. We have already shown excellent simulations of the mixing process (Figure 6). For this purpose, the extent of homogeneity was taken the same as those in the respective experimental measurements which varied over a wide range of 90 –99%. It was thought desirable to understand the relationship between the mixing time and the extent of homogeneity. For this purpose, the simulations have been performed for the case of Weiland and Onken (1981) and Bello et al. (1984). The results are given in Table 8. In all the cases the mixing time was found to be equivalent to the response time of a second order system (Figure 5).

Comparison Between Predicted and the Experimental Data for Axial Dispersion Coefficient Figure 5. A typical response curve at the outlet of riser for mixing time estimation.

To obtain the axial dispersion coefficient in the riser, the riser part, along with the connection of downcomer to the

Table 3. Comparison of CFD predictions with the experimental data of Weiland and Onken (1981).

VG (m s21)

VL (m s21)

0.006 0.014 0.020 0.035 0.038 0.046

0.11 0.16 0.21 0.25 0.28 0.32

(0.005) (0.012) (0.020) (0.030) (0.036) (0.045)

(0.128) (0.175) (0.222) (0.261) (0.294) (0.320)

G [ 0.007 0.016 0.027 0.033 0.041 0.049

(0.008) (0.017) (0.026) (0.034) (0.042) (0.050)

DL (m2 s21)

Material balance


Energy balance



LHS-1

RHS-1

LHS-2

RHS-2

Exp.

CFD

Exp.

CFD

nt (m2 s21)

0.60 0.69 0.77 0.80 0.85 0.89

0.55 0.65 0.72 0.73 0.79 0.84

0.023 0.028 0.035 0.040 0.042 0.046

0.013 0.023 0.03 0.040 0.039 0.046

0.005 0.009 0.012 0.014 0.016 0.017

0.0044 0.008 0.011 0.016 0.018 0.020

632.00 429.20 335.17 296.47 262.26 241.38

620 415 345 310 260 255

0.00014 0.00017 0.00020 0.00031 0.00043 0.00059

umix (sec)

The bracketed numbers indicate the experimental values.
LHS-1 ¼ V¯G/[ G; RHS-1 ¼ C0(V¯G þ V¯L) þ C1.

LHS-2 ¼ volume integral of 1, RHS-2 ¼ energy input.

Trans IChemE, Part A, Chemical Engineering Research and Design, 2006, 84(A8): 677 –690

684

ROY et al. Table 4. Comparison of CFD predictions with the experimental data of Bello et al. (1984). Material balance


AD/AR

VG (m s21)

VL (m s21)

[G

Energy balance



umix (s)

LHS-1

RHS-1

LHS-2

RHS-2

Exp.

CFD

0.11

0.016 0.023 0.035 0.057 0.087

(0.014) (0.022) (0.035) (0.056) (0.086)

0.070 0.082 0.095 0.110 0.131

(0.073) (0.085) (0.099) (0.116) (0.134)

0.036 0.053 0.076 0.113 0.154

(0.035) (0.052) (0.077) (0.112) (0.155)

0.40 0.42 0.45 0.50 0.55

0.43 0.43 0.46 0.45 0.55

0.010 0.025 0.037 0.060 0.085

0.009 0.023 0.039 0.061 0.082

52.30 46.00 42.44 38.29 37.82

51.0 45.5 41.0 40.0 39.5

0.25

0.015 0.022 0.036 0.055 0.088

(0.014) (0.022) (0.035) (0.056) (0.086)

0.130 0.155 0.179 0.208 0.242

(0.134) (0.156) (0.182) (0.213) (0.245)

0.030 0.042 0.064 0.090 0.125

(0.030) (0.043) (0.063) (0.092) (0.126)

0.46 0.51 0.55 0.61 0.68

0.50 0.52 0.57 0.61 0.69

0.014 0.029 0.044 0.069 0.097

0.011 0.027 0.042 0.070 0.095

54.77 45.54 44.58 43.10 43.00

52.5 44.0 43.5 43.0 41.7

0.44

0.014 0.023 0.037 0.058 0.088

(0.014) (0.022) (0.035) (0.056) (0.086)

0.200 0.220 0.260 0.320 0.350

(0.200) (0.237) (0.276) (0.323) (0.373)

0.025 0.038 0.055 0.075 0.110

(0.025) (0.037) (0.053) (0.076) (0.103)

0.56 0.60 0.66 0.73 0.83

0.56 0.61 0.65 0.77 0.81

0.013 0.032 0.059 0.073 0.109

0.012 0.030 0.061 0.070 0.107

56.50 53.66 52.15 50.09. 48.63

58.0 54.5 53.0 51.0 50.0

The bracketed numbers indicate the experimental values.
LHS-1 ¼ V¯/[G; RHS-1 ¼ C0(V¯G þ V¯L) þC1;

LHS-2 ¼ volume integral of 1, RHS-2 ¼ energy input.

Table 5. Comparison of CFD predictions with the experimental data of Verlaan et al. (1989).

VG (m s21) 0.015 0.038 0.064 0.082 0.110 0.140

(0.014) (0.035) (0.056) (0.078) (0.104) (0.139)

VL (m s21) 0.20 0.30 0.38 0.42 0.46 0.50

(0.20) (0.31) (0.38) (0.43) (0.46) (0.48)

[G 0.029 0.050 0.080 0.100 0.120 0.140

(0.026) (0.050) (0.070) (0.090) (0.109) (0.136)

DL (m2 s21)

Material balance


Energy balance



LHS-1

RHS-1

LHS-2

RHS-2

Exp.

CFD

Exp.

CFD

nt (m2 s21)

0.54 0.70 0.80 0.87 0.95 1.02

0.55 0.70 0.76 0.84 0.91 0.98

0.017 0.071 0.110 0.120 0.130 0.190

0.014 0.078 0.100 0.110 0.120 0.170

0.035 0.036 0.037 0.038 0.039 0.040

0.033 0.034 0.035 0.041 0.046 0.050

93 76 70 65 62 61

95.4 74.4 67.4 66.40 63.4 56.2

0.00057 0.00071 0.001 0.0011 0.0014 0.0015

umix (s)

The bracketed numbers indicate the experimental values.
LHS-1 ¼ V¯G/[G; RHS-1 ¼ C0(V¯G þ V¯L) þC1;

LHS-2 ¼ volume integral of 1, RHS-2 ¼ energy input.

Table 6. Comparison of CFD predictions with the experimental data of Kawase et al. (1994). Material balance
AD/AR

VG (m s21)

VL (m s21)

[G

Energy balance



umix (s)

LHS-1

RHS-1

LHS-2

RHS-2

Exp.

CFD

0.46

0.100 0.020 0.030 0.039 0.050

(0.009) (0.018) (0.026) (0.035) (0.045)

0.101 0.110 0.130 0.140 0.150

(0.102) (0.110) (0.120) (0.134) (0.161)

0.021 0.039 0.058 0.078 0.090

(0.020) (0.040) (0.056) (0.071) (0.083)

0.45 0.45 0.46 0.49 0.54

0.46 0.51 0.50 0.51 0.50

0.014 0.022 0.034 0.044 0.066

0.015 0.018 0.033 0.048 0.062

29.30 22.20 20.62 20.44 20.27

28.30 21.5 21.0 21.0 19.0

0.20

0.010 0.019 0.029 0.038 0.050

(0.009) (0.018) (0.026) (0.035) (0.045)

0.053 0.060 0.061 0.068 0.080

(0.053) (0.060) (0.063) (0.069) (0.082)

0.025 0.046 0.070 0.086 0.103

(0.023) (0.047) (0.065) (0.084) (0.099)

0.40 0.38 0.40 0.42 0.45

0.41 0.42 0.41 0.44 0.50

0.007 0.011 0.015 0.021 0.500

0.006 0.008 0.010 0.020 0.030

23.18 19.00 16.00 15.88 15.70

23.0 20.0 17.5 17.0 17.0

The bracketed numbers indicate the experimental values.
LHS-1 ¼ V¯G/[G; RHS-1 ¼ C0(V¯G þV¯L) þ C1;

LHS-2 ¼ volume integral of 1, RHS-2 ¼ energy input.

riser, was simulated separately. The information of inlet liquid velocity to the riser was taken from the complete simulation of EL-ALR, which was used to estimate the mixing time. Residence time distribution curves were obtained from the CFD simulation. Figure 7 shows the typical outlet response curves for superficial gas velocity of

5 mm s21 (Weiland and Onken, 1981). For this purpose, concentrations were estimated at the topmost axial location and at all the radial locations. The residence time distribution curves were used for the estimation of axial dispersion coefficient following the stepwise procedure suggested by Levenspiel (1995). The mean solute

Trans IChemE, Part A, Chemical Engineering Research and Design, 2006, 84(A8): 677– 690

CFD SIMULATION OF FLOW AND AXIAL DISPERSION

685

Table 7. Comparison of CFD predictions with the experimental data of Gavrilescu and Tudose (1997). Material balance
AD/AR

VG (m s21)

VL (m s21)

[G

Energy balance



DL (m2 s21)

umix (s)

LHS-1

RHS-1

LHS-2

RHS-2

Exp.

CFD

Exp.

CFD

nt (m2 s21)

0.122

0.010 0.020 0.041 0.069 0.089 0.109 0.139

(0.010) (0.020) (0.040) (0.060) (0.080) (0.100) (0.120)

0.050 0.063 0.083 0.089 0.095 0.102 0.105

(0.050) (0.064) (0.083) (0.089) (0.095) (0.100) (0.106)

0.027 0.048 0.087 0.131 0.160 0.191 0.231

(0.027) (0.048) (0.087) (0.134) (0.163) (0.189) (0.226)

0.370 0.395 0.436 0.477 0.509 0.534 0.575

0.365 0.385 0.424 0.474 0.499 0.504 0.536

0.054 0.099 0.188 0.270 0.315 0.380 0.480

0.052 0.097 0.172 0.240 0.312 0.320 0.420

0.0050 0.0063 0.0088 0.0117 0.0146 0.0145 0.0153

0.0042 0.0065 0.0082 0.0105 0.0138 0.0151 0.0162

741.62 555.93 384.77 304.52 263.65 237.84 227.17

748.2 550.4 376.2 310.8 270 242.6 220.2

0.000090 0.00011 0.00032 0.00048 0.00073 0.00091 0.00105

0.040

0.010 0.020 0.039 0.061 0.081 0.10 0.126

(0.010) (0.020) (0.040) (0.060) (0.080) (0.100) (0.120)

0.021 0.028 0.037 0.040 0.044 0.047 0.049

(0.022) (0.029) (0.039) (0.041) (0.044) (0.047) (0.050)

0.030 0.060 0.11 0.140 0.180 0.210 0.249

(0.029) (0.056) (0.101) (0.143) (0.180) (0.211) (0.249)

0.344 0.357 0.386 0.419 0.444 0.473 0.502

0.327 0.356 0.387 0.419 0.436 0.456 0.478

0.047 0.090 0.162 0.231 0.277 0.327 0.360

0.050 0.092 0.137 0.226 0.257 0.328 0.365

0.0069 0.0108 0.0137 0.0195 0.0223 0.0242 0.0250

0.0070 0.0110 0.013 0.020 0.021 0.022 0.026

650.790 459.040 321.200 253.070 212.120 183.320 169.670

645.6 450.2 326.2 250.0 216.4 180.0 172.2

0.0000580 0.0000866 0.000148 0.000319 0.000437 0.000520 0.000610

The bracketed numbers indicate the experimental values.
LHS-1 ¼ V¯G/[G; RHS-1 ¼ C0(V¯G þ V¯L) þ C1;

LHS-2 ¼ volume integral of 1, RHS-2 ¼ energy input.

concentration at any height z is given by Ð 2p Ð R

c u r dr d u C ¼ Ð0 2p Ð0 R 0 0 u r dr d u

(24)

When the reactor is modelled as a closed vessel, the first and second moments of the residence time distribution curve are related to the mean residence time (t), and the axial dispersion coefficient, DL, respectively (Levenspiel, 1995) as

Figure 6. Parity plot for mixing time: B, Weiland and Onken (1981); D, Bello et al. (1984); O,Verlaan et al. (1989); O, Kawase et al. (1994).

P tC(t)Dt tffi P C(t)Dt P 2 t C(t)Dt s2 ffi P  t2 C(t)Dt
 s2 DL DL 2 2 su ¼ 2 ¼ 2 (1  evL H=DL ) t vL H vL H vL H 2 s DL ffi 2 u

(25) (26) (27) (28)

Table 8. Effect of degree of homogeneity on mixing time.

umix (s) Authors

VG (m s21)

VL (m s21)

Exp.

Weiland and Onken (1981)

0.006 0.014 0.020 0.035 0.038 0.046

(0.005) (0.012) (0.020) (0.030) (0.036) (0.045)

0.11 0.16 0.21 0.25 0.28 0.32

(0.128) (0.175) (0.222) (0.261) (0.294) (0.320)

632.00 429.20 335.17 296.47 262.26 241.38

Bello et al. (1984) (AD/AR ¼ 0.44)

0.014 0.023 0.037 0.058 0.088

(0.014) (0.022) (0.035) (0.056) (0.086)

0.200 0.220 0.260 0.320 0.350

(0.200) (0.237) (0.276) (0.323) (0.373)

56.50 53.66 52.15 50.09 48.63

CFD (90%)

CFD (95%)

620 415 345 310 260 255 34.5 32.0 30.0 28.0 27.5

The bracketed numbers indicate the experimental values.

Trans IChemE, Part A, Chemical Engineering Research and Design, 2006, 84(A8): 677 –690

870 555 485 443 393 339 43 41 40.5 39 37

CFD (99%) 1000 685 610 580 543 465 58.0 54.5 53.0 51.0 50.0

686

ROY et al.

Figure 7. RTD curve of riser for VSG ¼ 5 mm s21 (Weiland and Onken, 1981).

where, t is the mean residence time, DL is the axial liquid dispersion coefficient, H is the height of the column and s represents the variance. The simulations have been carried out for experimental data of Weiland and Onken (1981), Verlaan et al. (1989), Gavrilescu and Tudose (1997) and Vial et al. (2005). The comparison between the predicted results and experimental data is given in Tables 3, 5, 7 and 9. Except Vial et al. (2005), all the axial dispersion values were predicted within 10% (Figure 8). In the case of Vial et al. (2005) the CFD under predicts the axial dispersion coefficient (Table 9). It can be seen from Tables 3, 5, 7 and 9 that there exists a wide variation of results in the published literature. In the

Figure 8. Parity plot for dispersion coefficient: S, Weiland and Onken (1981); þ, Verlaan et al. (1989); D, Gravrilescu and Tudose (1997) [AD/AR ¼ 0.12]; †, Gravrilescu and Tudose (1997) [AD/AR ¼ 0.04]; O, Vial et al. (2005) [AD/AR ¼ 0.25]; B, Vial et al. (2005) [AD/AR ¼ 0.28].

low VG range, Verlaan et al. (1989) have reported relatively high values of axial dispersion coefficient. This happens mainly due to the 908 connection of downcomer to the riser. This connection gives rise to the oscillation of gas plume through out the riser (Figure 9), which gives rise to recirculation like bubble column, resulting in high axial dispersion coefficient. The geometry of Gavrilescu and Tudose (1997) also has similar kind of connection, in spite of that they reported axial dispersion coefficient which is comparable with the other reported values. As shown in Figure 10, the plume oscillation in the riser dies

Table 9. Comparison of CFD predictions with the experimental data of Vial et al. (2005). Material balance
AD/AR

VG (m s21)

VL (m s21)

Energy balance



DL (m2 s21)

[G

LHS-1

RHS-1

LHS-2

RHS-2

Exp.

CFD

nt (m2 s21)

0.25

0.010 0.020 0.040 0.066 0.085 0.140 0.150 0.180 0.200 0.220 0.260

(0.009) (0.019) (0.037) (0.063) (0.084) (0.130) (0.150) (0.170) (0.190) (0.210) (0.250)

0.20 0.25 0.29 0.34 0.37 0.41 0.42 0.41 0.43 0.45 0.46

(0.20) (0.23) (0.28) (0.33) (0.36) (0.40) (0.41) (0.42) (0.43) (0.44) (0.47)

0.038 (0.040) 0.058 (0.056) 0.084 (0.082) 0.10 (0.11) 0.12 (0.13) 0.17 (0.16) 0.16 (0.17) 0.17 (0.18) 0.18 (0.19) 0.19 (0.20) 0.21 (0.22)

0.23 0.34 0.45 0.57 0.65 0.81 0.88 0.94 1.00 1.05 1.14

0.27 0.37 0.49 0.60 0.70 0.85 0.97 1.0 1.04 1.11 1.20

0.008 0.012 0.022 0.035 0.043 0.100 0.110 0.150 0.170 0.210 0.280

0.007 0.011 0.021 0.035 0.041 0.090 0.110 0.160 0.180 0.220 0.240

0.015 0.017 0.018 0.020 0.022 0.025 0.031 0.035 0.037 0.044 0.052

0.005 0.007 0.007 0.008 0.008 0.010 0.020 0.028 0.030 0.034 0.048

0.0000464 0.0000499 0.0000552 0.000144 0.000154 0.000438 0.000654 0.000865 0.000979 0.00105 0.00112

0.28

0.012 (0.011) 0.05 (0.023) 0.039 (0.037) 0.047 (0.047) 0.071 (0.071) 0.089 (0.088) 0.097 (0.096) 0.11 (0.113)

0.19 0.24 0.28 0.30 0.34 0.36 0.37 0.39

(0.17) (0.22) (0.26) (0.28) (0.33) (0.35) (0.36) (0.37)

0.035 (0.037) 0.061 (0.060) 0.080 (0.081) 0.093 (0.094) 0.11 (0.12) 0.13 (0.13) 0.14 (0.14) 0.15 (0.15)

0.30 0.38 0.46 0.50 0.59 0.68 0.69 0.75

0.32 0.40 0.48 0.51 0.62 0.70 0.72 0.75

0.028 0.031 0.045 0.062 0.13 0.12 0.13 0.15

0.030 0.032 0.050 0.067 0.11 0.13 0.14 0.15

0.022 0.023 0.025 0.022 0.031 0.038 0.051 0.061

0.013 0.015 0.017 0.018 0.017 0.024 0.028 0.039

0.0000328 0.0000596 0.0000642 0.000071 0.000157 0.000183 0.00019 0.00025

The bracketed numbers indicate the experimental values.
LHS-1 ¼ V¯G/[G; RHS-1 ¼ C0(V¯G þ V¯L) þC1;

LHS ¼ volume integral of energy dissipation, RHS ¼ energy input.

Trans IChemE, Part A, Chemical Engineering Research and Design, 2006, 84(A8): 677– 690

CFD SIMULATION OF FLOW AND AXIAL DISPERSION

Figure 9. Liquid velocity in EL-ALR, VG ¼ 56 mm s21 (Verlaan et al., 1989).

687

Figure 11. Liquid velocity in a part of EL-ALR, VG ¼ 63 mm s21, AD/ AR ¼ 0.25 (Vial et al., 2005).

down within a short distance from the bottom entry. This is mainly due to the fact that the inlet liquid velocity to the riser was considerably lowered by increasing the diameter of the connection as compared to the geometry of Verlaan et al. (1989). Furthermore, the height of the riser for the geometry of Gavrilescu and Tudose (1997) is higher than that of Verlaan et al. (1989), which is another reason for the decay of plume oscillation after a certain height. That height onwards the flow become gas –liquid co-current upflow which is characteristically close to pipe flow, which justifies the low value of axial dispersion in compare

to Verlaan et al. (1989). In contrast, as seen in Figures 11 and 12, for both the geometry of Vial et al. (2005) the flow is almost gas –liquid co-current up flow like in the case of Gavrilescu and Tudose (1997), but still the reported value of axial dispersion coefficient is high in comparison to other reported values. We have not understood the reasons of high value of DL for this case.

Figure 10. Liquid velocity in a part of EL-ALR,VG ¼ 60 mm s21, AD/AR ¼ 0.04 (Gavrilescu and Tudose 1997).

Figure 12. Liquid velocity in a part of EL-ALR, VG ¼ 71 mm s21, AD/AR ¼ 0.28 (Vial et al., 2005).

Trans IChemE, Part A, Chemical Engineering Research and Design, 2006, 84(A8): 677 –690

688

ROY et al.

Effect of AD/AR Ratio on Mixing Time and Axial Dispersion AD/AR ratio also has a prominent effect on mixing time (Figure 13), because the AD/AR ratio has a strong influence on the liquid velocity. A decrease in AD/AR decreases the liquid velocity in the riser. One can expect much higher mixing times for low AD/AR (lower liquid velocity) than at high AD/AR ratio. However, a decrease in the AD/AR ratio means a decrease in the total liquid volume which decreases the mixing time. Furthermore, for low AD/AR a larger fraction of the total liquid volume is contained in the riser section than for high AD/AR. Since riser section is the actively sparged zone an increase in the riser liquid volume to the total liquid volume ratio could promote more mixing and will cause lower overall mixing time. Effect of AD/AR ratio on axial dispersion is shown in Figure 14. The decreased AD/AR ratio increases the axial dispersion in riser due to the decrease in liquid velocity, which makes the flow in riser more like bubble column.

eddy diffusion to the mixing and axial dispersion coefficient. The effect of convection (i.e., liquid velocity) on mixing is shown in Figure 14. With an increase in VG, liquid velocity also increases and as a result mixing time decreases. In order to understand the effect of eddy diffusion on mixing time, numerical experiments were performed out, in which mixing time curve were simulated with different eddy diffusivity values. It can be seen from Figure 15 that the effect of eddy diffusivity on mixing time is negligible. Average eddy diffusivity values in riser are given in Tables 3, 4, 5, 7 and 9. It can be seen in all the cases that the values of axial dispersion coefficient is at least

Relative Contribution of Convection and Eddy Diffusion on Mixing Time and Axial Dispersion Coefficient In the published literature a confusion still persists (untill 2005) regarding the relative contribution of convection and

Figure 14. Effect of superficial gas velocity on dispersion coefficient with experimental data of Gavrilescu and Tudose (1997): S, AD/AR ¼ 0.04; O, AD/AR ¼ 0.12; —, CFD predictions.

Figure 13. Effect of AD/AR on mixing time (a) Kawase et al. (1994): S, AD/AR ¼ 0.46; B, AD/AR ¼ 0.20; —, CFD predictions. (b) Bello et al. (1984): S, AD/AR ¼ 0.11; O, AD/AR ¼ 0.25; B, AD/AR ¼ 0.44; —, CFD predictions.

Figure 15. Sensitivity of eddy diffusivity on mixing time, VG ¼ 14 mm s21 (Verlaan et al., 1989): (1) 10% (umix ¼ 90.4 s), (2) 50% (umix ¼ 93.6 s), (3) 100% (umix ¼ 95.4 s), (4) 150% (umix ¼ 97.4 s], (5) 200% (umix ¼ 97.6 s) (100% means CFD predicted eddy diffusivity).

Trans IChemE, Part A, Chemical Engineering Research and Design, 2006, 84(A8): 677– 690

CFD SIMULATION OF FLOW AND AXIAL DISPERSION

689

Table 10. Effect of liquid velocity on axial dispersion coeffiecient for VG ¼ 20 mm s21 (Weiland and Onken, 1981). DL (m2 s21) VL (m s21) 0.05 0.11 0.22 0.27 0.35

%

[G

nt (m2 s21)

Exp.

CFD

% change

22.7 50 100 122.72 159.1

0.060 0.055 0.050 0.047 0.040

0.00023 0.00019 0.002 0.00022 0.00027

0.012 0.012 0.012 0.012 0.012

0.084 0.038 0.011 0.0028 0.0024

600 216.6 28.3 276.66 280.0

100% means CFD predicted VL.

Table 11. Effect of eddy diffusivity on axial dispersion coeffiecient (VG ¼ 5 mm s21, VL ¼ 128 mm s21) (Weiland and Onken, 1981). DL (m2 s21) nt(%)

Exp.

CFD

10 50 100 150 200

0.005 0.005 0.005 0.005 0.005

0.0044 0.0044 0.0044 0.0044 0.0044

10 times higher than those of eddy diffusivity. This observation indicates eddy diffusivity cannot be a sole contributor to the axial dispersion coefficient. In order to understand the relative effect of convection and eddy diffusion on axial dispersion, two additional numerical experiments were performed, one with varying liquid velocity for same VG (Table 10) and the other with varying eddy diffusivity from 10% of the original value to 200% (Table 11). It can be seen from Table 10 that, as liquid velocity decreases for same VG the riser approaches bubble column behaviour of complete back mixing and gives 600% higher value of axial dispersion coefficient with only 22.7 times low value of liquid velocity, whereas the change in average eddy diffusivity is insignificant. When liquid velocity increases for same VG riser approaches plug flow behaviour and gives low value of axial dispersion coefficient. This observation indicates the overriding contribution of convection on axial dispersion, which again can be confirmed from the results presented in Table 11, which shows no effect of eddy diffusivity in axial dispersion.

CONCLUSIONS (1) The CFD model has been used to simulate the hydrodynamics in the riser of external airlift reactor and an extensive comparison of predicted and the experimental axial liquid velocity and hold-up profile has been presented. The agreement has been shown to be good. (2) A comparison of predicted and the experimental liquid circulation velocity has been presented over a wide range of VG, AD/AR and AR. The agreement has been shown to be excellent. (3) The CFD model has been successfully extended for prediction of mixing time for the experimental geometries covered by Weiland and Onken (1981), Bello

et al. (1984), Verlaan et al. (1989) Gavrilescu and Tudose (1997). Again a good agreement was observed. (4) The relative contribution of convection and eddy diffusion on mixing time and axial dispersion was studied. It was found that the axial dispersion is mainly controlled by the convection and the eddy diffusivity has a negligible role on the axial dispersion as well as mixing time. (5) The CFD model has been successfully extended for prediction of axial dispersion coefficient for the experimental cases of Weiland and Onken (1981), Verlaan et al. (1989), Gavrilescu and Tudose (1997), Vial et al. (2005). A good agreement was observed except Vial et al. (2005).

NOMENCLATURE A Bo c cˆ C C00 CB CD CL C0, C1 CTD Cm Cm,BI C11 C12 dB do D Deff DL Dm Dt EI g G H k MI MD ML MVM MTD N P r R Re t u Vo

area of cross section Bodenstein number (¼vLH/DL) instantaneous concentration of the tracer, kmol m23 Reynolds average concentration, kmol/m3 mean concentration of the tracer, defined in equation (24) initial concentration of the tracer interface energy transfer factor drag force coefficient lift force coefficient drift flux constants turbulent dispersion coefficient constant in k–1 model constant in bubble induced turbulence model model parameter in turbulent dissipation energy equation model parameter in turbulent dissipation energy equation bubble diameter, m sparger hole diameter, m diameter of the column, m effective diffusion coefficient, m2 s21 axial dispersion coefficient, m2 s21 molecular diffusivity, m2 s21 eddy diffusivity, m2 s21 rate of energy released by all the bubble in the column, J s21 gravitational constant, m s22 generation term height, m turbulent kinetic energy per unit mass, m2 s22 total interfacial force acting between two phases, N m23 drag force, N m23 lift force, N m23 added mass force acting, N m23 turbulent dispersion force, N m23 number of sparger hole pressure, N m22 radial distance, m column radius, m Reynolds number (2dBVS/n) time, s velocity vector, m s21 velocity of gas through sparger hole, m s21

Trans IChemE, Part A, Chemical Engineering Research and Design, 2006, 84(A8): 677 –690

690 V VS vL z

ROY et al. superficial velocity, m s21 axial slip velocity between gas and liquid, m2 s21 average liquid velocity, m s21 axial distance along the column, m

Greek symbols [ fractional phase holdup [ average fractional phase holdup 1 turbulent energy dissipation rate per unit mass, m2 s23 m molecular viscosity, Pa s21 meff effective viscosity, Pa s21 mBI bubble induced viscosity, Pa s21 mT turbulent viscosity, Pa s21 r density, kg m23 u tangential co-ordinate umix mixing time, s 2 s variance, defined in equation (12), s2 sk Prandtl number for turbulent kinetic energy s1 Prandtl number for turbulent energy dissipation rate s2u dimensionless variance, defined in equation (14) t mean residence time, defined in equation (11), s tk shear stress of phase k, Pa Subscripts D disp G k L R

downcomer dispersion gas phase phase, k ¼ G: gas phase, k ¼ L: liquid phase liquid phase riser

REFERENCES Akita, K., Okazaki, T. and Koyama, H., 1988, Gas holdup and friction factors of gas-liquid tow-phase flow in an air-lift bubble column, J Chem Eng Japan, 21: 476–482. Bello, R.A., Robinson, C.W. and Moo Young, M., 1984, Liquid circulation and mixing characteristics of airlift contactors, Can J Chem Eng, 62: 573– 577. Bentifraouine, C., Xuereb, C. and Riba, J.P., 1997, An experimental study of the hydrodynamic characteristics of external loop airlift contactors, J Chem Tech Biotech, 69: 345–349. Bendjaballah, N., Dhaouadi, M., Poncin, S., Midoux, N., Hornut, J.M. and Wild, G., 1999, Hydrodynamics and flow regimes in external loop airlift reactors, Chem Eng Sci, 54: 5211–5221. Chisti, M.Y. and Moo Young, M., 1988, Liquid circulation in ALR, Chem Eng Sci, 43: 451–457. Dhotre, M.T. and Joshi J.B., 2004, Two-dimensional CFD model for prediction of pressure drop and heat transfer coefficient in bubble column reactors, Trans IChemE, Part A, Chem Eng Res Des, 78: 689 –707. Deen, N.G., Solberg, T. and Hjertager, B.H., 2001, Large eddy simulation of gas– liquid flow in a square cross-sectioned bubble column, Chem Eng Sci, 56: 6341– 6349. Ekambara, K. and Joshi, J.B., 2003, CFD simulation of residence time distribution and mixing in bubble column, Can J Chem Eng, 81: 669 –676. Gavrilescu, M. and Tudose, R.Z., 1997, Mixing studies in external-loop reactors, Chem Eng J, 66: 97–104. Hunt, J.C.R., Auton, T.R., Sene, K., Thomas, N.H. and Kowe, R., 1987, ICHMT International Seminar on Transient Phenomena in Multiphase Flow, Dubrovnikm Yugoslavia, 103 –125. Hsu, Y.C. and Dudukovic, M.P., 1980, Gas holdup and liquid recirculation in gas-lift reactors, Chem Eng Sci, 35: 135–141.

Jackobsen, H.A., Sannaes, B.H., Grevskott, S. and Svendsen, H.F., 1997, Modelling of vertical bubble driven flow, Ind Eng Chem Res, 36: 4052–4074. Joshi, J.B., Ranade, V.V., Gharat, S.D. and Lele, S.S., 1990, Sparged lop reactors, Can J Chem Eng, 68: 705–741. Joshi, J.B., 2001, Computational flow modeling and design of bubble column reactors, Chem Eng Sci, 56: 5893–5993. Kawase, Y., Omori, N. and Tsujimura, M., 1994, Liquid-phase mixing in external-loop airlift bioreactor, J Chem Tech Biotech, 61: 49–55. Kochbeck, B., Lindert, M. and Hempel, D.C., 1992, Hydrodynamics and local parameter in three phase flow in airlift-loop reactors of different scale, Chem Eng Sci, 47: 3443–3450. Levenspiel, O., 1995, Chemical Reaction Engineering, 2nd edition (Wiley, New York). Lopez de Bertodano, M., 1991, Turbulent bubble flow in a triangular duct, PhD Thesis, Rensselaer Polytechnic Institute, Troy, New York. Merchuck, J.C. and Stein, Y., 1981, Local hold-up and liquid velocity in airlift reactors, AIChE J, 32: 1585–1596. Merchuck, J.C., 1986, Gas hold-up and liquid velocity in a twodimensional air lift reactor, Chem Eng Sci, 41: 11– 16. Nichol, R.S., 1984, Circulating bubble column, PhD Thesis, University of Cambridge, UK. Okada, K., Shibano, S. and Akagi, Y., 1993, Turbulent properties in bubble-flow region in external-loop airlift bubble column, J Chem Eng Japan, 6: 637–643. Popovic, M. and Robinson, C.W., 1987a, The specific interfacial area in external loop airlifts and a bubble column-I. aqueous sodium sulphate solution, Chem Eng Sci, 42: 2811–2824. Popovic, M. and Robinson, C.W., 1987b, The specific interfacial area in external loop airlifts and a bubble column-II. carboxymethyl cellulose/sulphite solution, Chem Eng Sci, 42: 2825–2832. Popovic, M. and Robinson, C.W., 1989, Mass transfer studies of external loop-airlifts and a bubble column, AIChE J 35: 393 –405. Sato, Y. and Sekoguchi, K., 1975, Liquid velocity distribution in twophase bubbly flow, Int J Multiphase Flow, 2: 79–95. Siegel, M.H., Merchuck, J.C. and Schugerk, K., 1986, Airlift reactor analysis: interrelationship between riser, downcomer and gas liquid separator, AIChE J, 27: 377 –388. Sokolichin, A., Eigenberger, G. and Lapin, A., 2004, Simulation of buoyancy driven bubbly flow: established simplifications and open questions, AIChE J, 50: 24 –45. Thakre, S.S. and Joshi, J.B., 1999, CFD simulation of bubble column reactors: Importance of drag force formulation, Chem Eng Sci, 54: 5055–5060. Verlaan, P., Tramper, J., Van’t Riet, K. and Luyben, K.Ch.A.M., 1986, A hydrodynamic model for an airlift loop bioreactor with external loop, Chem Eng J, 33: B43– B53. Verlaan, P., Tramper, J., Van’t Riet, K. and Luyben, K.Ch.A.M., 1989, Estimation of axial dispersion in individual sections of an airlift-loop reactor, Chem Eng Sci, 44: 1139– 1146. Vial, Ch., Poncin, S., Wild, G. and Midoux, N., 2005, Experimental and theoretical analysis of the axial dispersion in the liquid phase in external loop airlift reactors, Chem Eng Sci, 60: 5945–5954. Weiland, P. and Onken, U., 1981, Fluid dynamics and mass transfer in an airlift fermentor with external loop, Ger Chem Eng, 4: 42–50. Young, M.A., Carbonell, R.G. and Ollis, D.F., 1991, Airlift bioreactor: analysis of local two phase hydrodynamics, AIChE J, 37: 403–428. Zhang, D.Z. and Vanderheyden, W.B., 2002, The effects of mesoscale structures on the disperse two-phase flows and their closures for dilute suspensions, Int J Multiphase Flow, 28: 805–822. Zuber, N. and Findlay, J.A., 1969, Average volumetric concentration in two phase flow systems, J Heat Trans, 87: 453– 468. The manuscript was received 2 August 2005 and accepted for publication after revision 31 March 2006.

Trans IChemE, Part A, Chemical Engineering Research and Design, 2006, 84(A8): 677– 690