Experimental and Numerical Simulation for Gas− Liquid Phases Flow

Ind. Eng. Chem. Res. 2007, 46, 7317-7327

7317

Experimental and Numerical Simulation for Gas-Liquid Phases Flow Structure in an External-Loop Airlift Reactor Changqing Cao, Shuqin Dong, and Qingjie Guo* College of Chemical Engineering, Qingdao UniVersity of Science and Technology, Qingdao 266042, People’s Republic of China

A gas-liquid external-loop airlift reactor with a riser 0.47 m in diameter and 2.5 m in height and two externalloop down-comers 0.08 m in diameter and 2.5 m in height were used to investigate the gas-liquid two-phase flow structure. Local phase holdups were measured simultaneously by a microconductivity probe with air as the gas phase and water as the liquid phase over a wide range of operation conditions. Liquid flow velocity measurements were performed using the electrolyte tracer measurement (ETM) technique. The hydrodynamics near the sparger zone, riser disengagement zone (zone 1), junction zone (zone 2), and down-comer disengagement zone (zone 3) were systematically examined using the CFDs at the local scale and at the riser scale, respectively. The simulation results showed that zones 1, 2, and 3 exhibit three different flow regimes, which were the secondary mixed flow regime, the mixed flow regime, and the homogeneous bubble regime, respectively. It was also indicated that turbulent kinetic energy and turbulent kinetic energy dissipation rate were influenced by a gas sparger. These results were necessary to explain these different regimes using computational fluid dynamics (CFD) to provide deeper insight at the local scale for reactor geometry, such as gas sparger, junction and disengagement zones as well as the gas-liquid phase flow microstructure. The simulation results at the local scale were difficult to obtain by experiment. The numerical simulating results of local gas holdup and local gas and liquid velocities agreed well with the experimental data at a low gas flow rate. However, large errors occurred in the simulations at a high gas flow rate, because of poor estimation of the influence of bubble-induced turbulence or the higher density of the tracer and the poor mesh refinement. The flow structure and turbulence parameters of the phases presented here were useful for designing gasliquid external-loop airlift reactors. Introduction Airlift reactors are one of the important classes of modified bubble columns. They are characterized by their main advantages, which are their easy construction and low consumption, and mixing induced by gas aeration.1 However, airlift reactors also have additional advantages over bubble columns, because of the higher liquid overall velocity and the higher intensity of turbulence. The external-loop airlift reactor (EL-ALR) constitutes another important class of such modified bubble columns in which the overall circulation of a continuous phase is induced by a differential aeration between the riser and the down-comer. The external circulation of liquid is expected to promote the overall mixing (and back-mixing) in the air-lift reactor; thus, higher liquid velocities should be related with higher liquid circulation and backflow rates in the full reactors. However, higher liquid velocity reduces local liquid circulation and, therefore, back-mixing. This gives flatter liquid velocity profiles, reducing the shear, which improves life preservation of microorganisms in biological applications.2 Because the fractional gas holdup in the different sections is the critical parameter that determines the reactor hydrodynamics, most studies about airlift reactors had been concerned with estimating the variation of gas holdup in each section of the reactor with superficial gas velocity (UG).3-5 The superficial circulation velocity of the liquid phase is often identified as the key parameter for design. Consequently, many models and correlations had been proposed to predict the liquid velocity in airlift reactors.5-8 It seems that global parameters had been * To whom correspondence should be addressed. Tel.: 0086-53284022506. E-mail: [email protected].

widely studied but that little attention was paid to the local characteristics of the flow. Indeed, the influence of the local phenomena has great importance in regard to reactor performance, as shear stress, mixing, and mass transfer result simultaneously from the presence of an overall liquid circulation and from the local fluctuations with time of the hydrodynamic parameters. These are mainly due to bubble passage and bubble wake.9 Okada et al.10 and Young et al.11 proposed an experimental description of the radial evolution of the local gas holdup and of the local axial velocity of both phases. Little information on Reynolds radial normal and shear stress was available in the literature.12 For airlift loop reactors, experimental evidence is generally more recent, but less abundant, and it seems often contradictory (compare, e.g., the results of Gavrilescu and Tudose13 with those of Merchuk et al.14). This stems from the fact that the role of many geometrical parameters, such as the gas sparger, had often been ignored, as shown by Merchuk and Yunger15 for internalloop reactors and Snape et al.16 for EL-ALR. The influence of reactor geometry had also been studied. First, the effect of the cross-sectional area ratio of down-comer to riser had been shown to have a strong influence on the circulation velocity.5,17,18 The shape of the different section was also known to affect the performance of airlift reactors, especially in EL-ALR, which enables the use of wider range of geometry for each section.19 In contrast, most investigators have shown that the gas sparger, which was of important in conventional bubble columns,20,21 had little influence on the hydrodynamics in airlift reactors.22,23 No information about the effect of the junction zone and disengagement zone geometry on the local and global performance of EL-ALR was available in the literature.

10.1021/ie070690g CCC: $37.00 © 2007 American Chemical Society Published on Web 09/22/2007

7318

Ind. Eng. Chem. Res., Vol. 46, No. 22, 2007

Table 1. Summary of the Grids Used in the Computational Fluid Dynamics (CFD) Calculations grid type, GT

grid number, GN

orifice surface grid type, OSGT

orifice surface grid number, OSGN

non-orifice surface grid number, NOSGN

18

19440

Riser Sparger heterozygous tetrahedron

47374

wedge

18960

triangle Down-comer

The lack of data in the literature on the local flow characteristics is also true in the field of computational fluid dynamics (CFD). Generally, only classical one-dimensional (1-D) models were used up to now, based on the energy balance or the momentum balance equations. Both approaches give the same results, but the momentum balance method is often preferred.11,24,25 Merchuk and Stein26 and Verlaan et al.27 have also developed hydrodynamic models based on the drift flux model of Zuber and Findlay.28 Most of the models proposed for bubble columns are inadequate at high superficial liquid velocity. A few contributions studied airlift reactors using CFD: only a few works are available (Sokolichin and Eigenberger,29 Cockx et al.,30 Cockx et al.,31 and van Baten et al.32), but these authors simulated the behavior of an internal-loop airlift reactor (IL-ALR). Meanwhile, Vial et al.33 simulated the hydrodynamics in the riser of an EL-ALR using CFD. To obtain further information about the hydrodynamics in the riser of EL-ALR, the effect of airflow rate and reactor geometry such as sparger, junction zone, and disengagement zone on gas holdup, local liquid velocity, global turbulent kinetic energy, and turbulent kinetic energy dissipation rate is systematically studied. In this study, a microconductivity probe that measures technique and the electrolyte tracer measurement (ETM) technique are used to measure local phase holdups and liquid flow velocity with higher accuracy, respectively. The global hydrodynamics description includes different regimes, overall gas holdup in the riser, and liquid circulation velocity. Knowledge of global hydrodynamics also includes turbulent kinetic energy and turbulent kinetic energy dissipation rate values, which cannot yet be measured by simple methods, but can be evaluated using CFD. The local hydrodynamics description presented here consists of the radial evolution of local gas holdup and the time-averaged velocity of the gas and liquid phases. The experimental data are compared with numerical simulations, using a commercial CFD code (FLUENT 6.0) and classical assumptions of CFD models from the literature. Thus, it is shown that the current possibilities and limitations of CFD are able to predict both phenomena at the global and the local scales in EL-ALR and explain the discrepancies between the calculations and the experiments.

equation is an adopted renormalized group theory (RNG) k- model. Numerical Models. The continuity equations for continuous as well as dispersed phases in this case become

∂RqFq u qF q) ) 0 + ∇‚(Rqb ∂t where 2

∑ Rq ) 1 q)1 The momentum balance becomes, in general formulation,

F

(

)

[(

For the CFD calculations, the Eulerian-Eulerian approach has been chosen, because of its obvious computational advantages at high dispersed phase contents. Whereas, in EulerianLagrangian computations, the computational effort grows with increasing dispersed phase holdup. With gas holdup (up to 20 vol %), an Eulerian-Lagrangian calculation would yield high computational demands, which are unattainable on standard personal computer (PC) workstations. A three-dimensional reactor model is used as the physical model. The turbulent

)

∂ui ∂ui ∂ui ∂uk ∂P ∂ + Uk )+ µ + ∂t ∂xk ∂xi ∂xk ∂xk ∂xi

]

Fu′iu′k + FFi (2a) Ui ) ui + u′i

(2b)

P ) p + p′

(2c)

where Ui, ui, and u′i denote the instantaneous velocity, timeaveraged velocity, and fluctuating velocity, respectively. P, p, and p′ are the instantaneous pressure, time-averaged pressure, and fluctuating pressure, respectively. Liquid-phase turbulence has been used in the RNG k- model, as described by Yakhot et al.34 As a two-equation turbulence model, the RNG k- model introduces two additional variables into the calculations, namely, the turbulent kinetic energy (k) and the turbulent kinetic energy dissipation rate (). The RNG k- model is in the most general formulation, as follows. For the k equation:

F

(

)

∂k Dk ∂ R µ + Gk + Gb - F - YM + Sk (3) ) Dt ∂xi k eff ∂xi

where

Gk ) µkS2 Theoretical and Numerical Models

(1)

S ≡ x2Sij Sij Sij )

(

) ( )

∂uj ∂ui ∂uj + ∂xi ∂xj ∂xi

Gb ) -βgi

µt ∂T Prt ∂xi

YM ) 2FMt2

Ind. Eng. Chem. Res., Vol. 46, No. 22, 2007 7319

where R is the sound velocity.

Mt )

x

k R2

For the  equation:

F

(

)

D ∂  ∂ Rµ + C1 (Gk + C3Gb) ) Dt ∂xi  eff ∂xi k C2F

2 - R + S (4) k

where for Cµ ) 0.0845, η0 ≈ 4.38, β ) 0.012, and η ≡ Sk/.

µt ) R)

{

FCµk2 

}

CµF η3[1 - (η/η0)]2 2 k 1 + βη3

The effective liquid viscosity (µeff) in eqs 3 and 4 is then given by

d

( ) (x F 2k ) 1.72 xµ

νˆ

νˆ 3 - 1 + Cν

)

dνˆ

(5)

where Cν ≈ 100 and νˆ ) µeff/µ. The model parameters C1,C2, and C3 in eq 4 have been set constant: C1) 1.42, C2) 1.68, and C3 ) 1.30. Numerical Methods. A commercial CFD code (FLUENT 6.0) is used to simulate the local flow properties in EL-ALR. To solve the system of partial and ordinary differential equations previously presented numerically, discretion of the equations has been made using a finite volume scheme with an algebraic multigrid solver (AMG) as implemented in the CFD code CFX4.3. The coupling between pressure and velocity is addressed using the SIMPLEC method. The second discretion method is performed to solve each phase velocity weight, local phase holdup, and turbulent weight. The momentum equation coefficient has been implemented by means of the weighting method. Variables in the cell nodal point are simulated using second wind interpolation. Variables are stored using the center point method. Because of the large dimensions of the flow domain under consideration (a riser with an EL-ALR height of 2.5 m and diameter of 0.47 m; a down-comer with an EL-ALR height of 2.5 m and diameter of 0.08 m), the need is for transient and three-dimensional calculations, as well as limited computational resources. A very coarse numerical grid must be implemented with a total number of 66334 cells. The time step length is chosen as 0.5 s. A typical solver running over 10 s of computed time takes ∼100 h to complete on a supercomputer. Calculations always assume fully fluidized state as an initial condition. Convergence is good in all computations presented here, reaching the desired accuracy in 30 or less iterations. The ultimate convergent standard for all equations is no change of the volume fraction in the flow field. Boundary Conditions. In this study, an entire EL-ALR is taken as the computational domain. The divided network cells are carried out by primarily using a heterocomplex rectangular body structure in the riser and a wedge-shaped body structure in the downer. The gas distributor consists of 46 perforated upward orifices 0.002 m in diameter. The surface of the inlet gas orifice uses a triangle grid, and each surface of orifice is

Figure 1. (a) Schematic representation of the experimental airlift reactor. Legend: 1, manometer; 2, rotameter; 3, valve; 4, gas sparger; 5, conductivity meter; 6, A/D converter; 7, computer; 8, riser; 9, down-comer disengagement zone; 10, riser disengagement zone; 11, down-comer; 12, measuring electrodes; 13, tracer injection; 14, conductive probe; 15, junction zone. (b) Schematic diagram of the multiple-orifice sparger.

Figure 2. Schematic diagram of the three-wire electrode probe. Legend: 1, probe; 2, insulation film; 3, restle; 4, static wire; and 5, fill tube for tracer.

divided into 18 grids. To divide the grid conveniently, the outlet gas orifice is located at a distance of 0.01 m. The threedimensional (3-D) computational model mesh, the gas distributor with an orifice mesh, and the gas distributor with an orifice direction upward mesh are summarized in Table 1. The aforementioned meshes are yielded using GAMBIT software. Because of the presence of bend phenomenon intensively in EL-ALR, the wall surface is examined using a nonequilibrium wall function. Boundary conditions in the inlet have been adopted by entrance gas superficial velocity (the value is directly transferred on the basis of gas flux). Initial superficial gas velocities (UG0) from 0.02 m/s to 0.35 m/s are given for the entire simulation. The inlet turbulent kinetic energy (k0) is defined as

7320

Ind. Eng. Chem. Res., Vol. 46, No. 22, 2007

k0 ) 1.5(uI)2 where

u ) uz I ) 0.16(ReDt)-1/8 and the turbulent kinetic energy dispersed rate (0) is defined as

0 )

Cµ0.75k01.5 0.075Dt

The boundary conditions in the outlet have been used at the pressure outlet. Because the outlet position is a free interface, the meter pressure (Pm) is assumed to be zero.

Figure 3. Comparison of the calculated computational fluid dynamics (CFD) and the experimental technique for the estimation of VL radial profiles for various superficial gas velocities.

Experiments Experimental Setup. Figure 1 shows a schematic diagram of the EL-ALR used in this study. The column is composed of organic glass, which enables video recording and direct observation. The setup mainly consists of a external-loop upward pipe (0.47 m in diameter and 2.5 m in height), two externalloop downward pipes (0.08 m in diameter and 2.5 m in height), a gas distributor, a gas-liquid separating chamber, a gas phase hold-up measurement system, a liquid velocity measurement system, two dished valves, and a gas supply system. Gas is fed through a multiple-orifice distributor with 56 orifices 0.002 m diameter and with an upward orifice direction, uniformly spaced on the surface of the distributor. The superficial gas velocity (UG) is varied from 0.02 m/s to 0.35 m/s by means of a pair of rotameters. All experiments are performed with tap water at room temperature and ambient pressure. The liquid level in the riser disengagement zones is assumed to be sufficiently high for a complete deaeration, to avoid the presence of gas in the down-comer. A valve, located at the down-comer, is used to modify the decrease in pressure that is due to friction effects. The valve has an effect that is similar to reducing the diameter of the down-comer. Five valve openings (VOs)s100%, 80%, 60%, 40%, and 20%shave been used in this study. Measuring Instruments. The local time-averaged gas holdup measurement system includes a two-tip conductivity probe, a conductance meter (LDD Model 501), an amplifier, an A/D converter, and a computer. The analog signal is digitized by the data acquisition system at a frequency of 100 Hz. The local time-averaged gas holdup is then calculated by counting the fraction of the total measuring time when the tip is in contact with the gas phase. The structure of the two-tip conductivity probe is represented as seen in the work of Changqing et al.35 The circulation velocity of the liquid phase in the downcomer is determined using an invasive measurement technique called electrolyte tracer measurement (ETM). A pulse of 10 mL of saturated NaCl solution is injected at the top of the downcomer. The signals from three-wire electrode probes located in the down-comer are recorded simultaneously at a frequency of 200 Hz. Because the liquid phase exhibits almost-plug flow behavior, the conductivity signals can be fitted with a normal function by means of the least-squares method. From the difference between the moments of the normal distributions, the liquid velocity in the down-comer can be deduced. The superficial liquid velocity in the riser is then estimated

Figure 4. Comparison of calculated CFD and the experimental technique for the estimation of VG radial profiles for various superficial gas velocities.

Figure 5. Comparison of calculated CFD and the experimental technique for estimation of local gas holdup radial profiles for various superficial gas velocities.

using a sample mass balance on the liquid phase. The measurement principle is based on mass transfer at a probe surface influenced by the liquid velocity close to that surface. During measurements, a highly constant voltage is applied between the electrode surface (made from thin platinum wire) and a reference electrode made from thin copper plate. Increasing the liquid flow velocity yields a decreasing boundary-layer thickness at the electrode surface, which leads to increased mass transfer and subsequently increased electric current at constant voltage. After careful calibration, three-wire needle probes can deliver two-dimensional liquid flow velocity fields. Figure 2 shows the ETM three-wire needles electrode that has been used in this project. Several measurements have been performed at the same

Ind. Eng. Chem. Res., Vol. 46, No. 22, 2007 7321

Figure 6. Time-averaged velocity vector plot for the gas phase near the sparger. (Z ) 0.05-0.40 m, UG ) 0.234 m/s, VO ) 60%.)

Figure 8. Time-averaged velocity vector plot for the liquid phase near the sparger. (Z ) -0.15-0.15 m, UG ) 0.234 m/s, VO ) 60%.)

Figure 7. Time-averaged velocity vector plot for the gas phase near zone 3. (Z ) 2.22-2.50 m, UG ) 0.234 m/s, VO ) 60%.)

Table 2. Comparison of the CFD Calculated and Measured Global Hydrodynamics Parameters Simulations

Experiments

UGr (cm/s)

RGr

ULd (m/s)

RGr

ULd (m/s)

3.8 6.9 12.6 18.3

0.047 0.073 0.12 0.17

-1.43 -1.72 -2.05 -2.51

0.058 0.095 0.18 0.25

-1.4 -1.67 -2.00 -2.23

location, to check this point and to improve the accuracy of the time-averaged values. Results and Discussion Comparison of Simulation and Experimental Technique for Estimation Values. There is good agreement at the lower gas velocity between the experimental technique for estimation and the calculated global hydrodynamic parameters of superficial gas velocity (UGr) and superficial liquid velocity (ULr) in the riser. A comparison is reported in Table 2 when VO is 100%. As UGr increases, the difference between the calculations and experimental technique for estimations increases. Consequently, the calculations overestimate superficial liquid velocity (ULd) in the down-comer and underestimate the volume-

Figure 9. Time-averaged velocity vector plot for the liquid phase near zone 2. (Z ) 1.7-2.5 m, UG ) 0.234 m/s, VO ) 60%.)

averaged gas holdup (RGr) in the riser at a high gas flow rate. However, these results can be used even at high UGr values, as a rough estimation of the global hydrodynamics. In contrast, larger discrepancies appear at the local scale. This may be because the local liquid velocity (νLr) and local gas velocity (νGr) values predicted with the CFD model in the riser are expected to be different for a fine grid, whereas the high density of salt tracers is usually associated with significant deviations of residence and circulation time measurements. Respective comparisons between the experimental technique for estimations and simulated local liquid velocity (νLr) and local gas velocity (νGr) values in the riser are reported in Figures 3 and 4. These results show that the experimental technique for estimations and CFD computations agree only at lower gas flow rates, which is consistent with most of the works that have applied CFD to bubble columns until now. The discrepancies in Figures 3 and 4 are considered to be the effect of tracer density and the coarse mesh, as is probably the case for higher superficial gas velocities. Figure 5 shows that similar results are obtained with local gas holdup: the discrepancies increase as UGr increases. An open question is the influence of turbulence formulations on CFD accuracy, as a homogeneous single-phase turbulence

7322

Ind. Eng. Chem. Res., Vol. 46, No. 22, 2007

Figure 10. Cross-sectional time-averaged liquid velocity vector map with UG and axial position in the riser. (a) Z ) 0.8 m, UG ) 0.304 m/s; (b) Z ) 1.2 m, UG ) 0.342 m/s; (c) Z ) 1.4 m, UG ) 0.304 m/s; and (d) Z ) 2.0 m, UG ) 0.342 m/s.

is assumed and a strongly anisotropic bubble-induced turbulence has been shown to prevail. This is a critical point, because the pressure drop and, consequently, the overall liquid circulation is estimated from the k value. In the same way, local liquid velocity profiles are dependent directly on the estimation of k. As a consequence, a good prediction of turbulent parameters can be achieved only if k from the k- model and from the analytical model of Reynolds shear stress agrees. The results show that the k values from CFD computations are probably underestimated, because bubble-induced turbulence in the normal directions is not taken into account. As a consequence, the pressure drop is underestimated, which leads to overestimation of the overall liquid velocity and volumeaveraged gas holdup at higher gas flow rates. This conclusion is in agreement with Table 2. These results are also demonstrated that inadequate turbulent information is mainly responsible for the poor predictive ability of the CFD code at higher gas flow rates. Gas-Phase Hydrodynamics near the Sparger, Riser Disengagement Zone (Zone 1), Junction Zone (Zone 2), and Down-Comer Disengagement Zone (Zone 3). Figure 6 shows the result for the time-averaged gas-phase velocity vector plot near the sparger with UG ) 0.234 m/s and a VO of 60%. Although the inlet concentration of the gas phase is stable, the gas-phase velocity near the sparger is not stably distributed. The locally high concentrations represent the bubble swarms. The formation and disaggregation of bubble swarms are observed from CFD simulations. The spatio-temporal flow structure of bubble swarms is dependent on the local velo-

cities and concentration of the gas phase. It is shown that there are lower gas velocities in the center regions than near the wall above the sparger (Z < 0.3 m). At a height above the sparger (Z g 0.3 m), however, the gas velocities increase toward the wall, reach a maximum, and then decrease at the walls. The gas flows from the center to the wall and changes with time. These results show that there is probably a certain degree of liquid back-mixing and unstable and complex bubble swarms flow near the sparger region in the EL-ALR. This is in agreement with classical results obtained in airlift reactors.11,33 The results for the time-averaged gas-phase velocity vector plot near zone 3 with UG ) 0.234 m/s and a VO of 60% is shown in Figure 7. At the experimental process, the liquid level in zone 1 is chosen to be sufficiently high for a complete deaeration, to avoid the presence of gas in the down-comer. The simulation results in Figure 7 show that most of the gas overflows from zone 1. However, there is still a small quantity of gas that enters into zone 3. It can be also seen that the flow behaviors of zones 1, 3, and 2 are closer to a secondary mixedflow regime, a homogeneous flow regime, and a mixed flow regime, respectively. The three zones exhibit three different flow regimes, which can be explained by gas expansion due to significant friction and hydrostatic pressure effect. These simulation results are very important, because the liquid velocity and the geometry structure of the junction have significant influence on different regimes when a complete deaeration is not achieved. These results are not consistent with the previous investigation;33,36 by supposition, a complete deaeration was

Ind. Eng. Chem. Res., Vol. 46, No. 22, 2007 7323

Figure 11. Liquid-phase turbulent kinetic energy profiles for three simulation runs with UG and VO ) 80%: (a) UG ) 0.234 m/s, (b) UG ) 0.304 m/s, and (c) UG ) 0.342 m/s.

achieved at zone 1. Therefore, it seems necessary to explain these different regimes using CFD, to provide deeper insight at the local scale for reactor geometry, such as the sparger, junction zone, and disengagement zone. Liquid-Phase Hydrodynamics near the Sparger, Zone 1, Zone 2, and Zone 3. Figure 8 shows the result for the timeaveraged liquid-phase velocity vector plot near the sparger with UG ) 0.234 m/s and a VO of 60%. It can be seen that strong liquid-phase back-mixing arises near the sparger zone. As a consequence, the mixing liquid flows from the center to the wall and upward flows along the wall region, because of the interaction between the liquid flows from the down-comer bottom and downward flows of back-mixing liquid in the riser bottom. This may be due to such a high intensity of liquid fluctuations, which seem to be predominantly induced by bubbles at a low gas flow rate near the sparger zone. This means that, even at the highest superficial liquid velocity in the riser (ULr) studied, bubble-induced turbulence cannot be neglected. In agreement with the results obtained for sparger zone gas velocity fluctuation (Figure 6), the high intensity of liquid velocity fluctuations near the sparger zone are dependent both on global liquid velocity (which can be adjusted by modifying the down-comer valves), and on local bubble characteristics, as well as interaction between the bubble and the distributing

plate orifice, which sponges essentially on the physicochemical properties of both two phases and the sparger structure. The result for the time-averaged liquid-phase velocity vector plot near zone 2 has been described in Figure 9, using CFD for UG ) 0.234 m/s and a VO of 60%. It can be found that there are maximum liquid velocity and microcirculation patterns at zone 2. This may be due to the sparger structure and junction zone configuration result in the boundary-layer separation, and the liquid velocity directions intensively change. As a consequence, maximum energy loss is observed at this zone. A comparison with the results of Vial et al.33 shows that the volume-averaged gas holdup and ULr are dependent widely on the entire reaction geometry, because a strong influence of the valve has been demonstrated. However, relatively little attention has been paid to the effect of the gas sparger structure and junction zone configuration in airlift reactors, and these parameters have been varied only rarely in experimental studies.16,37 No experimental results with the existences of maximum liquid velocity and microcirculation patterns result in the gas sparger structure and junction zone configuration have been reported at all. It can be also found that the flow behaviors of zones 1, 3, and 2 are more similar to that of a homogeneous flow regime, a secondary mixed-flow regime, and an axial diffusion flow regime, respectively. The three zones exhibit three

7324

Ind. Eng. Chem. Res., Vol. 46, No. 22, 2007

is negligible. The local liquid velocity profiles are also affected by reactor structure. Down-comer gas holdup in ELALR is usually very small. This fact has already proved in the Results and Discussion section. Thus, the differences between the riser and down-comer gas holdups seem to increase with the increase in axial positions, resulting in a higher liquid velocity. Hydrodynamics Analysis at the Riser Scale. Global hydrodynamics in the riser of EL-ALR have been investigated in detail by Vial et al.,8,33 as a function of the VO. Turbulent kinetic energy (k) and turbulent kinetic energy dissipation rate () were rarely studied in the literature. In this work, k and  can be defined as follows:

1 k ) u′iu′i 2

( )

 ) µt

Figure 12. Gas and liquid phases turbulent kinetic energies profiles for two simulation runs with UG ) 0.234 m/s and VO ) 80%: (a) gas-phase turbulent kinetic energy, and (b) liquid-phase turbulent kinetic energy.

different liquid-phase flow regimes, which are different from the gas-phase flow regimes that correspond to the same three zones. Cross-Sectional Liquid Velocity Vector Maps with Various Axial Positions and Gas Velocities. A cross-sectional timeaveraged liquid velocity vector map has been obtained using CFD and is reported in Figure 10 for UG ) 0.304 m/s (Figures 10a and 10c) and 0.342 m/s (Figures 10b and 10d) and VO value of 80% at various heights of the riser, respectively. Figure 10 shows that local liquid velocity profiles at the same cross section seem to be significant nonuniform. The local liquid velocities are relatively lower in the center and higher near the wall zone. Meanwhile, local liquid flow directions are upward flows in the center zone and downward flows near the wall zone. This reason may be that local liquid velocity profiles show bubble-induced turbulence prevails near the wall, whereas single-phase turbulence prevails in the center of the riser. Such a high intensity of liquid velocity fluctuations is induced essentially by bubble passage and bubble wake.9 It can be also seen in Figure 10 that the increase of airflow rate leads to an increase of the liquid velocity fluctuations, especially for low axial positions (Z e 1.2 m). However, for higher values of the axial positions (Z g 1.4 m), the effect of airflow rate

∂u′i ∂u′i ∂xk ∂xk

where the expressions of u′i and µt are the pulsating velocity and turbulent viscosity, respectively. The terms k and  are used to denote the average pulsating kinetic energy for unit mass fluid flow and cause a decrease in turbulent kinetic energy. Predictions of the liquid-phase turbulent kinetic energy profiles (kL) are shown in Figure 11 at the UG values of 0.234, 0.304, and 0.342 m/s and a VO of 80%. This figure shows that kL increases with UG, and its profiles present relative symmetry to the fluctuant mode with increasing UG. The result shows that liquid flow is unstable in EL-ALR, as shown in Figure 8. From Figure 11, it can also be observed that kL is comparatively influenced by the gas sparger. The value of kL is higher near the gas sparger zone. Apart from the gas sparger, kL is gradually reduced. This seems to prove that the occurrence of pulsating liquid is mainly related to local phenomena at the gas sparger zone. The profiles of gas-phase turbulent kinetic energy (kG) and kL in the riser are reported in Figure 12 at UG ) 0.234 m/s and a VO of 80%. The profiles of kG and kL seem to have almost the same contours and magnitudes. This may be due to intensive interactions between the gas phase and the liquid phase, as well as the strong anisotropy of turbulence, and thus it resulted in the stress tensor order of pulsating magnitude of kG and kL being at the same level. It is difficult to obtain a precise estimation of kG and kL, because these parameters are dependent both on the accuracy of pulsating velocity values of gas and liquid from measurements and on the validity of the model for kinetic viscosity estimation. Further work is needed to determine whether experimental measurement techniques used in this work allow concluding the deviations can be related with the poor mesh refinement due to the higher density of the tracer or with the underestimate bubble-induced turbulence with this particular CFD model may have a role. The profiles of the liquid-phase turbulent kinetic energy dissipation rate (L) in the riser are shown in Figure 13 at the respective UG values of 0.234, 0.304, and 0.342 m/s and a VO of 80%. It can be found that L is strongly dependent on UG. L increases as UG increases. This is consistent with the common vision that airlift reactors are characterized by a lower shear stress level, as compared to bubble columns. It may be expected that the turbulent intensity increases, as more momentum is transferred to the liquid phase, which could explain the increase in L when UG increases. Another reason may be that more mechanical energy is transformed to heat energy with increasing UG, and thus resulted in an increase in L.

Ind. Eng. Chem. Res., Vol. 46, No. 22, 2007 7325

Figure 13. Liquid-phase turbulent kinetic energy dissipation rate profiles for three simulation runs with UG and VO ) 80%: (a) UG ) 0.234 m/s, (b) UG ) 0.304 m/s, and (c) UG ) 0.342 m/s.

As shown in Figure 13, L is higher near the sparger zone. This means that the maximum gas velocities arise in the gas sparger outlet, and more mechanical energy is reduced by the friction process and thus results in bad turbulent kinetic energy dissipation. Thus, to obtain a savings in energy, we should consider preferentially design parameters of the gas sparger, because of mainly energy loss at the gas sparger outlet. Conclusions (1) A detail description of the local hydrodynamic parameters of both the liquid and gas phases has been obtained in an external-loop airlift reactor (EL-ALR), using the microconductivity probe measuring technique and the electrolyte tracer measurement technique. This may be useful for the purpose of computational fluid dynamics (CFD) validation. (2) All the experimental data have compared with simulations obtained using a classical CFD approach. A reasonable agreement is achieved in the homogeneous regime. The inaccuracy of simulations at high gas flow rates is shown to be due to the poor estimation of the turbulent parameters.

(3) The influence of reactor geometry on the local hydrodynamics has been simulated using the CFD code FLUENT 6.0. It is shown that there are lower gas velocities in the center region than near the wall above the sparger (Z < 0.3 m). At heights far above the sparger (Z g 0.3 m), however, the gas velocities increase toward the wall, reach a maximum, and then decrease at the wall. The mixing liquid flows from the center to the wall and upward flows along the wall region due to strong liquidphase back-mixing near the sparger zone. The simulation results also show that zones 1, 2, and 3 exhibit three different flow regimes, which are a secondary mixed-flow regime, a mixed flow regime, and a homogeneous bubble regime, respectively. Liquid velocity and the geometry structure of the junctions have significant influence on the different regimes. Meanwhile, maximum liquid velocity and microcirculation patterns exist at zone 2. Local liquid velocity radial profiles seem to be significantly nonuniform. The values are lower in the center zone than that near the wall zone. Local liquid flow directions are upward flows in the center zone and downward flows near the wall zone.

7326

Ind. Eng. Chem. Res., Vol. 46, No. 22, 2007

(4) Global hydrodynamics of k and  in the riser have also been simulated. It is shown that k increases and its profiles present from relative symmetry to fluctuant mode with increasing superficial gas velocities. The profiles of kG and kL seem to have almost the same contours and magnitudes; kG and kL are comparatively influenced by the gas sparger and are higher near the gas sparger zone. Apart from the gas sparger, kG and kL are gradually reduced. The value of G is strongly dependent on UG. The value of G increases as UG increases, and G is higher near the gas sparger zone.

F ) density (kg/m3) µ ) viscosity (Pa.s) µeff ) effective viscosity (Pa s) µt ) turbulent viscosity (Pa s) AbbreViations CFD ) computational fluid dynamics EL-ALR ) external-loop airlift reactor ETM ) electrolyte tracer measurement IL-ALR ) internal-loop airlift reactor VO ) valve opening

Acknowledgment The financial support of this research by Taishan Mountain Scholar Constructive Engineering Foundation of China (No. Js200510036), by National Natural Science Foundation of China (No. 20676064), and by Young Scientist Awarding Foundation of Shandong Province (No. 2006BS08002) is gratefully acknowledged. Nomenclature C1 ) constant in k- model in eq 4 C2) constant in k- model in eq 4 C3) constant in k- model in eq 4 Cµ ) constant in k- model in eq 4 Cν ) constant in k- model in eq 5 dp ) particle diameter (m) Dt ) column diameter (m) Gb ) turbulent kinetic energy produced terms by average velocity grads (kg/(m s3)) Gk ) turbulent kinetic energy produced terms by buoyancy (kg/ (m s3)) k ) turbulent kinetic energy (m2/s2) k0 ≡ initial turbulent kinetic energy (m2/s2) P ) pressure (Pa) P ) instantaneous pressure (Pa) p′ ) fluctuating pressure (Pa) Pm ) meter pressure (Pa) Re ) Reynolds number (dimensionless) Sk ) turbulent kinetic energy added term (kg/(m s3)) S ) turbulent kinetic energy dissipation rate added term (kg/ (m s4)) t ) time (s) UG ) superficial gas velocity (m/s) UGr ) superficial gas velocity in the riser (m/s) UL ) superficial liquid velocity (m/s) ULd ) superficial liquid velocity in the down-comer (m/s) ULr ) superficial liquid velocity in the riser (m/s) Ui ) instantaneous velocity (m/s) ui ) time-averaged velocity (m/s) u′i ) pulsating velocity (m/s) uz ) superficial velocity at z axial coordinate (m/s) VG ) local gas-phase velocity (m/s) VL ) local liquid-phase velocity (m/s) YM ) turbulent kinetic energy dissipation term produced by volume expansion (kg/(m s3)) Z ) axial position (m) R ) sound velocity (m/s) Rk ) constant in k- model (turbulent Prandtl number for k) R ) constant in k- model (turbulent Prandtl number for ) Rq ) qth phase holdup (dimensionless) β ) constant in k- model in eq 4 η0 ) constant in k- model in eq 4  ) turbulent kinetic energy dissipation rate (m2/s3)

Literature Cited (1) Deckwer, W.-D. Bubble Column Reactors; Wiley: Chichester, U.K., 1992. (2) Chisti, M. Y.; Moo-Young, M. Airlift Reactors: Characteristics, Applications and Design Considerations. Chem. Eng. Commun. 1987, 60, 195-242. (3) Siegel, M. H.; Merchuk, J. C.; Schugerl, K. Airlift Reactor Analysis: Interrelationships Between Riser, Down-Comer and Gas-Liquid Behavior, Including Gas Recirculation Effects. AIChE J. 1986, 32, 15851595. (4) Joshi, J. B.; Ranade, V. V.; Gharat, S. D.; Lele, S. S. Sparged Loop Reactors. Can. J. Chem. Eng. 1990, 68, 705-741. (5) Dhaouadi, H. Etude d’un re´acteur a` gazosiphon a` recirculation externe. Ph.D. Thesis, INPL, Nancy, France, 1997. (6) Popovic, C. W.; Robinson, C. W. External-Loop Airlift Bioreactors: Study of the Liquid Circulating Velocity in Highly Viscous NonNewtonian Liquids. Biotechnol. Bioeng. 1988, 32, 301-312. (7) Gavrilescu, M.; Tudose, R. Z. Study of the Liquid Velocity in External-Loop Airlift Reactor. Bioprocess Eng. 1995, 14, 183-193. (8) Vial, C.; Poncin, S.; Wild, G.; Midoux, N. Experimental and Theoretical Analysis of Axial Dispersion in the Liquid Phase in ExternalLoop Airlift Reactors. Chem. Eng. Sci. 2005, 60, 5945-5954. (9) Schmidt, J.; Nassar, R.; Lu¨bbert, A. Local Dispersion in the Liquid Phase of Gas-Liquid Reactors. Chem. Eng. Sci. 1992, 47, 203-213. (10) Okada, K.; Shibano, S.; Akagi, Y. Turbulent Properties in BubbleFlow Region in External-Loop Airlift Reactor. J. Chem. Eng. Jpn. 1993, 43, 637-643. (11) Young, M. A.; Carbonell, R.; Ollis, D. F. Airlift Bioreactors: Analysis of Local Two-Phase Flow Hydrodynamics. AIChE J. 1991, 37, 403-428. (12) Menzel, T.; In der Weide, T.; Staudacher, O.; Wein, O.; Onken, U. Reynolds Shear Stress for Modeling of Bubble Column Reactors. Ind. Eng. Chem. Res. 1990, 29, 988-994. (13) Gavrilescu, M.; Tudose, R. Z. Mixing Studies in External Loop Airlift Reactors. Chem. Eng. J. 1997, 66, 97-104. (14) Merchuk, J. C.; Contreras, A.; Garcia, F.; Molina, E. Studies of Mixing in a Concentric Airlift Reactor. Chem. Eng. Sci. 1998, 53, 709719. (15) Merchuk, J. C.; Yunger, R. The Role of Gas-Liquid Separator of Airlift in the Mixing Process. Chem. Eng. Sci. 1990, 45, 2973-2975. (16) Snape, J. B.; Zahradnik, J.; Fiavola´, M.; Thomas, N. H. Liquidphase Properties and Sparger Design Effect in an External Loop Airlift Reactor. Chem. Eng. Sci. 1995, 50, 3175-3188. (17) Bello, R. A. A Characterization Study of Airlift Contactors for Applications to Fermentations, Ph.D. Thesis, University of Science, Ontario, Canada, 1981. (18) Choi, K. H.; Lee, W. K. Circulation Liquid Velocity, Gas Holdup, Volumetric Oxygen Transfer Coefficient in External-Loop Airlift Reactors. J. Chem. Technol. Biotechnol. 1993, 56, 51-58. (19) Siegel, M. H.; Robinson, C. W. Applications of Airlift GasLiquid-Solid Reactors in Biotechnology. Chem. Eng. Sci. 1992, 47, 32153229. (20) Zahradnik, J.; Fiavola´, M.; Ruzicka, M.; Drahosˇ, J.; Kasˇta´nek, F.; Thomas, N. H. The Effect of Electrolytes on Bubble Coalescence and Gas Hold-up in Bubble Column Reactors. Chem. Eng. Sci. 1997, 52, 38113826. (21) Camarasa, S.; Vial, C.; Poncin, S.; Wild, G.; Midoux, N.; Bouillard, J. Influence of Coalescence Behavior of the Liquid and of Gas Sparging on Hydrodynamics and Bubble Characteristics in a Bubble Column. Chem. Eng. Process. 1999, 38, 329-344. (22) Bovosumbut, S.; Wilhelm, A. M.; Riba, J. P. Influence of Gas Distributor Design on the Oxygen Transfer Characteristics in an ExternalLoop Airlift Reactor. J. Chem. Technol. Biotechnol. 1987, 40, 167-174.

Ind. Eng. Chem. Res., Vol. 46, No. 22, 2007 7327 (23) Vial, C.; Camarasa, E.; Poncin, S.; Wild, G.; Midoux, N.; Bouillard, J. Study of the Hydrodynamic Behavior in Bubble Columns and External Loop Airlift Reactors Through Analysis of Pressure Fluctuations. Chem. Eng. Sci. 2000, 55, 2957-2973. (24) Glennon, B.; ai-Masry, W.; MacLoughlin, P. F.; Malone, D. M. Hydrodynamic Modeling in an Air-Lift Reactor. Chem. Eng. Commun. 1993, 121, 181-192. (25) Sa´ez, A. E.; Marquez, M. A.; Roberts, G. W.; Carbonell, R. G. Hydrodynamic Model for Gas-Lift Reactors. AIChE J. 1998, 44, 14131423. (26) Merchuk, J. C.; Stein, Y. A Distributed Parameter Model for an Airlift Fermentor. Effects of Pressure. Biotechnol. Bioeng. 1981, 13, 13091324. (27) Verlaan, P.; Tramper, J.; van’t Riet, K. A Hydrodynamic Model for an External Loop Airlift Reactor. Chem. Eng. J. 1986, 33, B43-B53. (28) Zuber, N.; Findlay, J. A. Average Volumetric Concentration in TwoPhase Flow Systems. Int. J. Heat Mass Transfer 1965, 87, 453-468. (29) Sokolichin, A.; Eigenberger, G. Gas-Liquid Flow in Bubble Columns and Loop Reactors: Part I. Detailed Modeling and Numerical Simulation. Chem. Eng. Sci. 1994, 49, 5747-5762. (30) Cockx, A.; Line´, A.; Roustan, M.; Do-Quang, Z.; Lazanova, V. Numerical Simulation and Physical Modeling of the Hydrodynamics in an Internal Airlift Loop Reactor. Chem. Eng. Sci. 1997, 52, 3787-3793. (31) Cockx, A.; Do-Quang, Z.; Line´, A.; Roustan, M. Use of the CFD for Simulating Hydrodynamics and Mass Transfer in Industrial Ozonation Towers. Chem. Eng. Sci. 1999, 54, 5091-5100.

(32) van Baten, J. M.; Ellenberger, J.; Krishna, R. Hydrodynamics of Internal Air-Lift Reactors: Experiments Versus CFD Simulations. Chem. Eng. Process. 2003, 42, 733-742. (33) Vial, C.; Poncin, S.; Wild, G.; Midoux, N. Experimental and Theoretical Analysis of the Hydrodynamics in the Riser of an ExternalLoop Airlift Reactor. Chem. Eng. Sci. 2002, 57, 4745-4762. (34) Yakhot, V.; Orszag, S.A. Renormaliazation Group Analysis of Turbulence I: Basic Theory. J. Sci. Comput. 1986, 1, 1-51. (35) Changqing, C.; Mingyan, L.; Qingjie, G. Experimental Investigation into the Radial Distribution of Local Phase Holdups in a Gas-LiquidSolid Fluidized Bed. Ind. Eng. Chem. Res. 2007, 46, 3841-3848. (36) Joshi, J. B. Computational Flow Modeling and Design of Bubble Column Reactors. Chem. Eng. Sci. 2001, 56, 5893-5933. (37) McManamey, W. J.; Wase, D. A. J.; Raymahasay, S.; Thayanithy, K. The Influence of Gas Inlet Design on Gas Holdup Values for Water and Various Solutions in Loop-type Airlift Fermenter. J. Chem. Technol. Biotechnol. 1984, B34, 151-164.

ReceiVed for reView May 15, 2007 ReVised manuscript receiVed July 24, 2007 Accepted August 6, 2007 IE070690G