Computational Fluid Dynamics Modeling and Validation of Bidisperse

6316

Ind. Eng. Chem. Res. 2008, 47, 6316–6326

Computational Fluid Dynamics Modeling and Validation of Bidisperse Fluidized Industrial Powders Olumuyiwa Owoyemi and Paola Lettieri* Department of Chemical Engineering, UniVersity College London, Torrington Place, London WC1E 7JE, U.K.

This work presents the validation of a computational fluid dynamics (CFD) model recently developed by Owoyemi et al. [Owoyemi, O.; Mazzei, L.; Lettieri, P. CFD modeling of bidisperse fluidized suspensions and investigation of the role of particle-particle drag on mixing and segregation. AIChE J. 2007, 53 (8), 1924-1940] for the simulation of bidisperse fluidized suspensions. The paper focuses on the fluidization of industrial rutile powders used in the titanium refining industry. Results are presented on the experimental and numerical investigation into the mixing and segregation behavior of a bidisperse fluidized mixture of natural rutile (flotsam) and slag (jetsam), which differ in size, having mean diameters of 186 and 305 µm, respectively, and having the same density of 4200 kg/m3. Three different average compositions, corresponding to the average mass fraction of jetsam particles of 0.25, 0.50, and 0.75 in the bed were considered. Computational results obtained for mixing and segregation patterns, bed expansion, and bubble dynamics are validated against experimental data. 1. Introduction Fluidized-bed technology is employed in a wide range of industrial applications, covering the pharmaceutical, food, chemical, and petrochemical industries as well as the mining and power generation industries. In most industrial applications of fluidization, the suspension consists of nonspherical particles of different diameters and sometimes different densities. Since the 1980s, when it was first realized that computational fluid dynamics (CFD) could be used for the quantitative design calculations for single-phase systems, its use has grown increasingly as a modeling tool for studying multiphase systems including fluidization. In process industries, for instance, the use of CFD in single-phase systems has led to reductions in the cost of product and process development activities (through the reduction of downtime), a reduction in the need for physical experimentation, shortened time to market, improved design reliability, and increased conversions and yields.1,2 In this work, the Eulerian-Eulerian modeling approach is used, which treats the particle and fluid phase as interpenetrating continua. In the modeling of monocomponent dense gas solid flows, the granular kinetic theory is typically employed by many researchers to translate the behavior of the many particles into one continuum.3,4 This approximation has, however, been proven to be inadequate in giving a fully theoretical description of the behavior of particles in dense gas-solid flows.5,6 This has led to the use of concepts borrowed from the field of soil mechanics as well as analogies drawn from the continuum theory of solid elasticity in providing a description adequate for the dense particulate phase.7–10 The picture gets even more complex when the EulerianEulerian approach is applied to modeling of dense bidisperse gas-solid systems. The problem of accounting for the behavior of the solid phases becomes more intricate. Constitutive expressions for the interaction between the two solid phases as well as suitable expressions for the behavior of each solid phase have to be provided. Some researchers have circumvented this issue by developing ingenious ways of reducing the complexity of the problem. One such way is modeling of the three-phase * To whom correspondence should be addressed. Tel: +44 20 7679 7867. Fax: +44 20 7383 2348. E-mail: [email protected].

system within the framework of a two-phase system.11 This reduces the problem to a two-fluid problem; however, care has to be taken to use suitable expressions to account for the drag force exerted by the fluid on the solid mixture as well expressions for the behavior of the mixed particles. Reviews on this subject are reported by Kuipers and van Swaaji,12 van den Akker,13 Sundaresan,14 Arastoopour,15 and Campbell.16 This work aimed at investigating the computational modeling of bidisperse gas-solid fluidized bed is concerned with the titanium refining industry where a bubbling fluidized bed is used for the extraction of titanium from naturally occurring ore. The refining process begins in a fluidized bed with the chlorination of titanium-rich rutile ore, which is composed of many constituents. Because of the size and density differences of all of the feedstock components used in the process, there are industrial concerns about the pervasiveness of dead zones within the fluidized bed as a result of feed stock segregation. Thus, the objective of the work is to develop a model capable of predicting the degree of mixing and segregation in the fluidizedbed system. To this end, the following powders, slag, and natural and synthetic rutile, belonging to the Geldart group B classification, were provided by Huntsman-Tioxide for the experimental and computational investigations in this project. A first step in the investigation showing the ability of the standard Eulerian granular kinetic theory approach available in the commercial code CFX to reproduce the fluidization behavior of the monocomponent natural and synthetic rutile powders was presented by Owoyemi et al.17 in this journal. The study on the monocomponent materials was further developed by Owoyemi and Lettieri,18 who presented a new expression for the fluid-particle interaction term, based on the fluidized-bed elasticity concept originally proposed by Wallis,19 together with a stability analysis and CFD validation. Owoyemi and Lettieri18 presented also a comparison between the new model and the granular kinetic theory model used in Owoyemi et al.17 and showed that collisional stresses had little impact on the numerical predictions for the materials investigated. The modeling approach developed by Owoyemi and Lettieri18 for monocomponent powders was subsequently extended to bidisperse systems using the three-fluid approach, where a

10.1021/ie800138n CCC: $40.75  2008 American Chemical Society Published on Web 06/25/2008

Ind. Eng. Chem. Res., Vol. 47, No. 16, 2008 6317

Figure 1. Particle size distribution of (a) natural rutile and (b) slag. Table 1. Continuity and Momentum Equations of Motion for a Bidisperse System continuity equation, fluid phase ∂ε + ∇ · (ε〈u〉f) ) 0 ∂t

(1)

∂φk + ∇ · (φk〈u〉Fk) ) 0 ∂t

(2)

continuity equation, solid phase k

momentum equation, fluid phase Ff

[ ∂t∂ (ε〈u〉 ) + ∇ · (ε〈u〉 〈u〉 ) ) ∇ · T - n 〈f 〉 f

f

f

f

f F1

1

- n2〈ff 〉F2 + εFfg

(3)

momentum equation, solid phase 1 F1

[ ∂t∂ (φ 〈u〉 ) + ∇ · (φ〈u〉 〈u〉 ) ) ∇ · T + n 〈f 〉

+ n1〈fp 〉F1 + φ1F1g

(4)

F2

[ ∂t∂ (φ 〈u〉 ) + ∇ · (φ〈u〉 〈u〉 ) ) ∇ · T + n 〈f 〉

- n1〈fp 〉F1 + φ2F2g

(5)

F1

F1

F1

1

1

f F1

1

momentum equation, solid phase 2 F2

F2

F2

2

separate momentum equation is solved for the fluid and each solid phase (see Owoyemi et al.1). The latter also analyzed the effect of the particle-particle drag force on the dynamics of the bidisperse system, by comparing three different closures available in the literature, namely, Gidaspow et al.,20 Syamlal,21 and Bell.22 The contribution of these closures was assessed against a reference test case where the particle-particle drag was not accounted for. Similar to the approach followed for the monocomponent systems, Owoyemi et al.1 presented a sensitivity analysis on the effect of the collisional stress on the simulations of the bidisperse systems and demonstrated that it did not affect the numerical prediction of the mixing index and coefficient of segregation. Owoyemi et al.1 validated their proposed model for a bidisperse system of two well-mixed monosized ballotini materials and showed the ability of the model to simulate mixing patterns for these materials. In this present paper, the model

2

f F2

2

presented by Owoyemi et al.1 is validated for bidisperse mixtures of industrial powders, i.e., natural rutile and slag powders, characterized by a size distribution as reported in Figure 1 and having mean diameters of 186 and 305 µm, respectively. The rutile powders investigated here also have a much higher density, equal to 4200 kg/m3, than the ballotini materials used in the previous work. This allows one to investigate the validity of the model for a wider range of materials’ characteristics. Furthermore, three different average compositions, corresponding to the average mass fractions of jetsam particles in the bed of 0.25, 0.50, and 0.75, were considered, and different fluidizing velocities were investigated in order to enable a thorough testing of the model for different degrees of mixing and segregation. The results provide qualitative and quantitative measures of the ability of the model to reproduce the behavior of bidisperse industrial materials. Computational results obtained for mixing

6318 Ind. Eng. Chem. Res., Vol. 47, No. 16, 2008

and segregation patterns, bed height, and bubble dynamics are validated against experimental data.

Table 2. Closure Equations for a Bidisperse System

2. CFD Model for Bidisperse Mixtures

fluid phase

The Eulerian-Eulerian averaged equations of motion used in this work are summarized, in tensorial form, in Table 1. The full details regarding the derivation of the equations of change and closure equations can be found in Owoyemi et al.1 The closure equations are reported in Table 2. Assuming a Newtonian behavior for each continuum, the effective stress tensors take the form as that in eqs 6 and 7 in Table 2. The internal stress associated with the solid phases, on the other hand, has not been modeled, as discussed in detail in Owoyemi et al.1 The fluid-particle interaction force accounts for the resultant force that is exerted by the fluid on the particles of phase Fk, where k refers to phases 1 and 2. In this model, this force is assumed to be made up of three contributions, namely, buoyant force, drag force, and elastic force. The buoyant force (see eq 8) associated with the phase Fk is related to the isotropic contribution to the effective fluid stress tensor and the solid volume fraction φk. The drag force is expressed, as usual, as the product of a drag coefficient βk and the relative slip velocity between the fluid and the particle phase of interest. Here the closure used for the interphase momentum transfer coefficient is that derived by Di Felice23 for the pressure drop per unit length in a suspension modified by Owoyemi and Lettieri18 (see eqs 9 and 11). The concept of elasticity associated with multiphase flow, and in particular with solid fluidized suspensions, was initially introduced by Wallis19 and thereafter further developed by Foscolo and Gibilaro9 and Gibilaro.24 A revised, multidimensional, and somewhat generalized formulation of the force has been subsequently advanced by Mazzei et al.25 for monocomponent systems. The force is therein regarded as a component of the drag force; more specifically, the contribution related to voidage gradients within the suspension that arise under nonequilibrium conditions. The approach of Mazzei et al.25 has been extended for bidisperse mixtures and presented in detail in Owoyemi et al.;1 the closure takes the expression as shown by eqs 12 and 14. The use of separate momentum equations for each particulate phase in the modeling of bidisperse mixtures requires an extra particle-particle drag force to account for the collisions between particles that belong to different size classes. Several expressions have been put forward for this.21,22,26–28 In general, the force is expressed as the product of a drag coefficient ζ and the relative slip velocity between the two solid phases (see eq 15 in Table 2). In this work, the correlation proposed by Syamal21 for the momentum transfer coefficient ζ is used (see eq 16 in Table 2). 3. Experimental Section The experimental setup used in this work consists of a 2D Plexiglas rectangular column, 600 mm high, 350 mm wide, and 10 mm thick. The distributor is a uniformly permeable sintered bronze rectangular plate with a thickness of 3.5 mm. The fluidizing gas, air, is supplied via rotameters. Pressure taps are installed 100 mm apart along the height of the bed from which pressure readings are collected via an electronic manometer. A system of two interlocked on/off valves operated simultaneously is installed on the rig to allow for instantaneous evacuation of the fluidizing gas during the bed freeze tests performed for the analysis of the mixing and segregation that occurs in the bed. Bed freeze analysis was conducted 1 min after the initiation of the fluidization experiments because this was the time indicated

Effective Stress Tensor

2 Tf ) -〈p〉fI + µf[∇〈u〉f + (∇〈u〉f)T] + (kf - µf)∇ · 〈u〉fI 3

(6)

solid phase 2 Tk ) -〈p〉FkI + µk[∇〈u〉Fk + (∇〈u〉Fk)T] + (kk - µk)∇ · 〈u〉FkI 3

(7)

Buoyancy Force nk〈f b 〉Fk ) -φk∇〈p〉f

(8)

nk〈f d 〉Fk ) βk(〈u〉f - 〈u〉Fk)

(9)

Drag Force

with βk )

[( ) 17.3 Rek

R

+ 0.336R

]

f Fk 1⁄R Ff|〈u〉 - 〈u〉 |φk -1.8

ε

dk

(10)

where Rek )

Ffε|〈u〉f - 〈u〉Fk|dk µf

and R ) 2.55 - 2.1[tanh(20ε - 8)0.33]3

(11)

Elastic Force nk〈f e 〉Fk ) Ek(∇ε · 〈nd 〉Fk)〈nd 〉Fk

(12)

where the elastic modulus Ek is

) 〉 | - φ (F - F )g

2 -3.8 + X |nk〈fd Ek ) dk 3 ε

[(

{

Fk

k

k

f

(13)

The functional form of X is

X)

[( ) ( ) ( ) ( )

1 1 - ln R R 17.3 Rek

17.3 Rek

R

ln

R

]

+ (0.336)R +

}

17.3 + (0.336)R ln(0.336) Rek {6.3[tanh(20ε 17.3 R + (0.336)R Rek ] }{sech (20ε - 8)

0.33 2

8)

2

0.33

}{6.6(20ε - 8)0.67}

(14)

Particle-Particle Interaction Force n1〈f p 〉F12 ) ς(〈u〉F2 - 〈u〉F1) ) -n2〈fp 〉F21

(15)

Syamlal21

(

3(1 + e) ς)

where

)

π π2 + Cf φ F φ F (d + d2)2go 2 8 1 1 2 2 1 |〈u〉F1 - 〈u〉F2| 2π(F1d13 + F2d23)

(

1 3 d1d2 φ1 φ2 go ) + 2 + ε ε d1 + d2 d1 d2

)

(16)

Ind. Eng. Chem. Res., Vol. 47, No. 16, 2008 6319

by Huntsman that was required to achieve steady-state conditions in their industrial fluidized-bed system. This information was used to perform the bed freeze test accordingly, in order to evaluate the degree of mixing achieved in the laboratory-scale fluidized bed after the same length of time. In the bed freeze tests, the fluidizing air supply was shut off abruptly. The bed at rest was then split into four equal horizontal layers, of about 7 cm in height in order to allow for a representative sample of the bed to be collected, and each layer was sieved to obtain the percentage by weight of the different components. Samples used for sieving in this work were collected by means of a probe attached to a vacuum pump. A schematic of the experimental apparatus can be found in Owoyemi et al.1 The materials used in this study are natural rutile and slag, which differ only in size, having mean diameters of 186 and 305 µm, respectively, and the same density of 4200 kg/m3. The mean particle diameter and size distribution were determined by means of sieving analysis. For each powder, the sieve test was repeated three times and average values for the mean particle diameter were taken, where the surface-volume mean particle diameter has been used to calculate the mean particle diameter. The particle size distribution for the powders investigated is shown in Figure 1. The natural rutile and slag particles represent the flotsam and jetsam particles, respectively. In all of the experiments performed, the bed is initially completely segregated. The flotsam particles are filled first and the jetsam particles are then added on top up to a height of 300 mm. The amounts of flotsam and jetsam particles that have been used for all experiments are summarized in Table 3. The mixture fractions of the powders were initially filled in two layers in which the flotsam particles (natural rutile, NR, smaller particles) occupied the bottom half of the bed and the jetsam particles (bigger particles, slag, SG) occupied the top half of the bed. The mixtures were then fluidized at a constant excess gas velocity, uxs. Here, the excess gas velocity is given by the following: uxs ) uo - umf, where uo refers to the superficial gas velocity and umf is the minimum fluidization velocity. The value of uxs was determined using the expression developed by Wu and Baeyens29 that relates the mixing index M representative of the bed mixedness to the visible bubble flow rate, as shown below:

( )

M ) 1 - 0.0067dR × 1.33

QB A

-0.75

(17)

Here dR is the ratio d1/d2 of the larger diameter to the smaller diameter of the particles. The visible bubble flow rate, QB/A, is given by the two-phase theory of fluidization by Davidson and Harrison:30 QB ) γ(uo - uinf) A

(18)

where the empirical correlation for γ advanced by Wu and Baeyens29 has been used: γ(Ar) )

{

1 for Ar e 50 -0.21 for Ar > 50 2.27Ar

(19)

The Archimedes number Ar of the mixture is given by Ar )

Fm - Ff Ga; Ff

1 ω1 ω2 ) + Fm F1 F2

and the Galileo number Ga is equal to

(20)

Ga)

dm3Ff2g 2

µf

; dm)

ω1F2 + ω2F1 dd ω1F2d2 + ω2F1d1 1 2

(21)

In eqs 20 and 21, the quantities ω1 and ω2 denote the mass fractions of the larger and smaller particles, respectively, in the powder devoid of fluid; these, therefore, take the expressions ω1)

φ1F1 φ2F2 ; ω2) φ1F1 + φ2F2 φ1F1 + φ2F2

(22)

The minimum fluidization velocity Umf of the mixture was calculated using the correlation proposed by Gossens et al.:31 Umf)

µf [(33.72 + 0.0408Ar)0.5 - 33.7] Ffdm

(23)

In summary, the excess gas velocity uxs was obtained using the following procedure: (1) a value for the mixing index M was chosen (equal to 0.9 in this work), and the visible bubble flow rate was evaluated from eq 17; (2) the minimum fluidization velocity umf of the mixture was then determined using eq 18; (3) the superficial gas velocity uo was finally assessed using eq 20. The superficial gas velocities employed for the various mixture fractions for both the experiments and simulations are summarized in the form of a “matrix” in Table 4. Here, the superficial gas velocities used for the investigations on mixing, for instance, can be discerned by selecting a particular column (jetsam mixture fraction) and thereafter selecting the desired excess gas velocity, in this case 0.20, and last adding the desired excess gas velocity to the minimum fluidization velocity specified in the chosen column. Digital video recordings of the fluidized bed were also made to analyze the development of the bubble dynamics within the bed and to determine the bubble size at the operating conditions employed. Images captured by means of a web camera at 14 frames/s, for 80 s, were recorded and subsequently analyzed using Optimas 6.0, an image analysis software. 4. Simulations In the present study, all of the simulations were carried out using a commercial CFD package, CFX 4.4. This package allows the free implementation of extra equations, boundary conditions, and differencing schemes. The closure relations reported in Table 2 were implemented in the code. As was previously pointed out, the internal stress associated with the solid phases has been neglected. The omission of this stress might give rise to localized overcompaction in some regions of the simulated fluidized bed. To counter this effect, the numerical algorithm excess solid volume correction developed by Owoyemi et al.1 for bidisperse systems has also been implemented within the CFX4.4 code. A thin 3D computational grid, 600 mm height × 350 mm width × 10 mm thickness, matching the size of the experimental setup used in this work, was used for the simulations. The front and back wall effects were neglected; this is consistent with the physical characteristics of the materials investigated, which did not display any interparticle force effects at the wall that could be visually observed. The left and right walls of the domain were modeled using no-slip velocity boundary conditions (for simplicity) for all phases. Dirichlet boundary conditions were employed at the bottom of the bed to specify a uniform gas inlet velocity. A pressure boundary condition was specified at the top of the bed and set to a reference value of 1.015 × 105 Pa. The distributor was made impenetrable to the

6320 Ind. Eng. Chem. Res., Vol. 47, No. 16, 2008 Table 3. Amounts of Flotsam and Jetsam Particles Used in the Experimental Study of the Industrial Powders ωjet

mjet (kg)

mflot (kg)

0.25 0.50 0.75

1.660 1.110 0.577

0.631 1.085 1.804

Table 4. Superficial Gas Velocities Used in the Study of Mixing and Segregation Patterns for Industrial Materials average mass fraction of jetsam, ωjet

quantity, m/s

0.25 uxs umf uo

0.50

0.20

0.10 0.060 0.260 0.160

0.75

0.20

0.10 0.078 0.278 0.178

0.20

0.10 0.107 0.307 0.207

Table 5. Computational Parameters Used in the CFD Simulations description gas density gas viscosity bed height settled bed height grid cell size time step superficial gas velocity coefficient of restitution coefficient of friction

symbol

value

Ff 1.29 µf 1.85 × 10-5 Hb 0.60 Hs 0.30 ∆x and ∆y 0.005 t 10-4 Uo 0.25 e 0.97 Cf 0.15

units

comments

kg/m3 Pa · s m m m square cells s m/s Syamlal21 Syamlal21

Table 6. Bidisperse Mixture Fractions Investigated simulation

ωflot

ωjet

1 2 3

0.25 0.50 0.75

0.75 0.50 0.25

solid phase. A second-order discretization scheme, SUPERBEE, was used for all equations in order to improve numerical convergence and the computational prediction of bubble shape and behavior (see Copper and Coronella32). This implicitly introduces a numerical dissipation term (or artificial viscosity), which is purely numerical in origin and which serves to improve and guarantee the stability of the solution (see Anderson33). A grid size of 0.5 cm and a time step of 1.0 × 10-4 s were adopted in this work. These values were chosen as a result of a sensitivity analysis carried out by Owoyemi et al.1 on the time step and grid resolution, where four simulations were carried out using grid sizes of 5 and 10 mm square cells as well as time steps of 10-4 and 10-3 s. The results were compared on the basis of macroscopic average fluidization properties such as the bed height, mixing index, and segregation for a bidisperse mixture of ballottini powders and showed that a grid size of 0.5 cm and a time step of 1.0 × 10-4 s were adequate and sufficiently small to match the experiments. The fluidization conditions used for all simulations are summarized in Table 5. Three different average compositions, corresponding to the average mass fraction of jetsam particles (slag) in the bed of 0.25, 0.50, and 0.75, were considered, as reported in Table 6. For each of the three cases reported in Table 6, the investigation was carried out at two excess gas velocities, uxs equal to 0.20 and 0.10 m/s, and each simulation was carried out for 10 s in real time. In all, a total of six numerical simulations were performed for this study. The overall computational time was roughly 85 days for each simulation. The simulations were carried out using three Dell Xeon P4 3.2 GHz dual processor machines.

Table 7. Experimental and Theoretical Values for Mixing Index M Obtained for Different Jetsam Fractions at an Excess Gas Velocity of 0.20 m/s mixing index M ωjet

CFD

experimental

Wu and Baeyens29

0.75 0.50 0.25

1.00 0.98 0.95

0.92 1.00 1.00

0.90 0.90 0.90

5. Results and Discussion 5.1. Mixing. The degree of mixing of the bidisperse mixture can be determined using the bed mixing index as defined by Rowe et al.34 as M)

〈ωjet 〉t 〈ωjet 〉o

(24)

where 〈ωjet〉t is the average mass fraction of the jetsam phase in the top region of the bed and 〈ωjet〉o is the average mass fraction of the jetsam phase evaluated over the entire bed. The top region of the bed is not uniquely defined but can be chosen somewhat arbitrarily. Here, following van Wachem et al.,11 the region was assumed to be the top 25% of the bed. The above definition intrinsically equates a state of perfect jetsam segregation at the bottom of the bed to a mixing index of M ) 0 and a state of perfect mixing to a mixing index of M ) 1. Perfect jetsam segregation at the top of the bed, which is usually imposed before the system is fluidized, is instead given by a mixing index of M ) (mjet + mflot)/mjet, where mjet and mflot are the overall masses of jetsam and flotsam in the bed, respectively, whose values for the different mixtures investigated are reported in Table 5. Table 7 reports a comparison of the mixing index M obtained from experiments with numerical results obtained after 10 s of simulation, using an excess gas velocity of 0.20 m/s, for all of the different jetsam fractions. The values of M are also compared with the prediction obtained from the semiempirical correlation developed by Wu and Baeyens.40 All of the simulations show a good agreement with the experimental data and the semitheoretical expression. In quantitative terms, results from the numerical simulations are within 8% of the experimental results and show a slight increase of the mixing index obtained from the experiments with an increase in the percentage of the jetsam compared to a slight decrease of the mixing index obtained from the CFD simulations. This discrepancy may be due to the pervasiveness of the mixing and solids recirculation mechanisms over the antagonistic effect of segregation in the CFD simulations, highlighting present gaps in current computer models to describe correctly the particle-particle interactions and, in turn, incorrectly predicting the mixing and segregation dynamics of these systems. A quantitative comparison between CFD and the prediction from the semiempirical correlation developed by Wu and Baeyens29 is found to be within 10%. However, the correlation seems to underpredict the experimental and computational mixing indexes in all of the cases investigated. Furthermore, as expected, the correlation gives a constant value for all three cases because it does not take into account the different average composition of the mixture but it is simply expressed as a function of the ratio between the larger and smaller diameters of the particles and the excess gas velocity, both of which are the same in all three cases. Nonetheless, a fluidized bed with a mixing index greater than 0.9 is usually considered to be in a well mixed state. Figure 2 shows an experimental snapshot of

Ind. Eng. Chem. Res., Vol. 47, No. 16, 2008 6321

Figure 2. Experimental snapshots of the fluidized bed using a mixture fraction of 0.75 slag at an excess gas velocity of 0.20 m/s.

Figure 3. Computational snapshots of the simulated fluidized bed using a mixture fraction of 0.75 slag at an excess gas velocity of 0.20 m/s.

the typical evolution of mixing in the fluidized bed for a mixture fraction, ωjet, of 0.75 jetsam (slag) at an excess gas velocity of 0.20 m/s. From the snapshot, it can be observed that the flotsam particles (natural rutile), which originally occupy the bottom part of the bed, are transported upward by rising bubbles. This mechanism of transport has also been identified by other researchers.35,36 In the present case, experimental observations showed that particles that gather in the wake of a rising bubble are carried from near the bottom of the bed to the surface of the bed. En route, natural rutile particles from the wake of the rising bubble are constantly exchanged with new particles from the dense surrounding region. This mode of suffusion of the natural rutile particles to the surface of the bed appears as the “fingerlike” propagation of the particles in the snapshot. At the same time, slag particles were seen to “rain down” through approaching bubbles as a mode of transport toward the bottom of the bed. The large presence of bubbles of different sizes in the fluidized bed adds complexity to the mixing pattern observed. Figure 3 shows a computational snapshot of the mixing mechanism in the simulated fluidized bed for a mixture fraction of 0.75 slag particles at an excess gas velocity of 0.20 m/s. From the figure, it can be seen that the experimentally observed mixing mechanism has been captured by the numerical simulation; natural rutile particles are observed to permeate to the top

Figure 4. Computational evolution of the jetsam mass fraction profile with time for a mixture concentration of 75% slag.

of the bed in a “fingerlike” manner in concurrence with earlier experimental findings. This transportation of natural rutile particles to the top of the bed is also facilitated by bubble motion, as captured by the simulation showing substantial agitation throughout the computational bed. The slag particles are also observed to travel downward at the walls of the vessel as well as in the wake of rising bubbles until evenly distributed.

6322 Ind. Eng. Chem. Res., Vol. 47, No. 16, 2008

Figure 5. Computational evolution of the jetsam mass fraction profile with time for a mixture fraction of 50% slag.

Figure 6. Computational evolution of the jetsam mass fraction profile with time for a mixture fraction of 25% slag.

Similar results were obtained for the other mixture compositions. Qualitatively, the experimental and computational results showed that a change in the average mixture of the jetsam fraction had no effect on the overall mixing mechanism in the bed. For all of the cases investigated, computational beds when fully mixed remained well mixed throughout. Figure 4 shows the evolution of the jetsam mass fraction profile with time obtained from the CFD simulations using a bidisperse mixture fraction of 0.75 slag particles. During the first 2 s, the jetsam phase spreads rapidly throughout the system; its mass fraction is seen to increase at the bottom of the bed until an even distribution is attained throughout the bed after about 4 s of simulation. There is little variation in the jetsam mixing profile thereafter, and therefore subsequent profiles are not reported. Figures 5 and 6 show the computational evolution of the slag particles with time for 50% and 25% jetsam concentrations. In both cases, the bed becomes fully mixed around 4 s. The total computational mixing time remains unchanged even with a change in the jetsam mixture fraction. This result may not be surprising because the excess gas velocity

was held constant at 0.20 m/s for all of the simulations, meaning that the same amount of gas is available for propagation into bubbles in all of the cases. 5.2. Macroscopic Bed Properties. A quantitative comparison of the bulk properties of the bed is shown in Table 8, in which results for the bed height, bed voidage, and bubble holdup are reported for simulations carried out at an excess gas velocity of 0.20 m/s. The properties have been time-averaged after the first 2 s of simulation in order to reduce the effect of perturbations associated with the bed startup. From Table 8, varying differences are observed between numerical predictions and experimental findings for the different jetsam mass fractions investigated. Quantitatively, a percent error spanning the range from 2.0% (25% jetsam fraction) to 4.0% (75% jetsam fraction) was found for bed height when computational and experimental results were compared. The slight deviation of computational predictions from experimental results might be due to the deficiency in accounting for a wide particle size distribution for the jetsam (slag) in the numerical model. A range between 5.6% (25% jetsam fraction) and 7.0% (75% jetsam fraction) error can be observed for the bed voidage in Table 8. Results show that the agreement between numerical and experimental results is better at lower jetsam mass fractions. A possible reason for this might be the modeling assumptions used in this work, where particles have been assumed to be spherical and, as a result, a bed voidage of 0.40 is intrinsically assumed in all of the simulations. The above assumption is in disparity with experimental observations, where the bed voidage at rest was higher. A comparison between numerical predictions and experimental results yielded a percent error in the range from 6.7% (25% jetsam fraction) to 97.8% (75% jetsam fraction) for the bubble holdup. The bubble holdup in this work is defined as the fraction of gas, above that required for minimum fluidization that leaves the fluidized bed through bubbles. Experimentally from Table 8, the bubble holdup was seen to increase with an increase in the slag mass fraction. Thus, it seems that an increase in the mass fraction of large particles causes excess gas to preferentially leave the bed via bubbles. Computationally, numerical values obtained for bubble holdup are constant and seem to only show good agreement at lower jetsam fractions, within 10%. This agreement becomes poorer with an increase in the jetsam concentration. This seems to be caused by the fact that computational bubbles obtained during this study were found to be smaller in size than the experimentally observed ones. A possible explanation for this may be found in a deeper understanding of the role played by the frictional viscosity used in the solid stress tensor on the mechanism of bubble growth in the fluidized-bed dynamics, hence requiring further investigation. On the other hand, it is also worth noting that the distributor plate may have an effect on the simulated bubble dynamics. In this work, the distributor plate is simulated as a fully open plate; using an orifice distributor instead may have an effect on both the size of the bubbles that form at the bottom of the bed and their growth.

Table 8. Comparison of Experimental and Numerical Predictions of Averaged Fluidized-Bed Macroscopic Quantities at an Excess Gas Velocity of 0.20 m/s bed height

bed voidage

bubble holdup

slag, ωjet

exp

CFD

% diff

exp

CFD

% diff

exp

CFD

% diff

0.75 0.50 0.25

37 37.5 38.5

35.5 36.2 37.7

4.0 3.5 2.0

0.537 0.545 0.552

0.499 0.510 0.521

7.0 6.4 5.6

0.360 0.236 0.190

0.182 0.175 0.178

97.8 34.9 6.7

Ind. Eng. Chem. Res., Vol. 47, No. 16, 2008 6323

Figure 7. Comparison of the simulated bubble diameter with the Darton et al.37 equation at an excess gas velocity of 0.20 m/s using a mixture fraction of 0.75 slag particles.

Figure 8. Comparison of the experimental bubble diameter with the Darton et al.37 equation at an excess gas velocity of 0.20m/s usi ng a mixture fraction of 0.75 slag particles.

5.3. Bubble Properties. The analysis of the bubble dynamics within the bed was carried out by comparing results obtained from numerical computations with data obtained from experimental investigations and with the Darton et al.37 bubble growth equation: db ) 0.54(u - umf)0.4(h + 4√A0)0.8 ⁄ g0.2

(25)

where db is the bubble diameter, h is the height above the distributor, A0 is the “catchment area”, which characterizes the gas distributor, and u is the superficial gas velocity. A value of 0 has been suggested by Darton et al.37 for A0 in the absence of available data on the distributor characteristics. A voidage of 0.80 was assumed to represent a computational bubble. The experimental analysis of the bubble diameter was conducted using Optimas 6.0. The computational bubble analysis, on the other hand, was performed using the numerical algorithm recently developed by Mazzei and Lettieri.38 Figure 7 shows a plot of the bubble diameter versus the bed height for all of the simulated bubbles computed at a constant excess gas velocity of 0.20 m/s using a jetsam concentration of 75%. The figure also shows predictions obtained from eq 19. Over 5000 bubbles, i.e., the total number of bubbles that has evolved during 10 s of real-time simulation, have been plotted on the graph. Figure 7 shows an enormous spread in predicted bubble diameters due to coalescence, breakup, and interaction of bubbles with the wall of the computational vessel. It can be seen that the predicted bubble diameters are smaller in the higher part of the fluidized bed and are within (20% of the predictions from the Darton et al.37 equation. The prominence of small bubbles in the bed may be due to less pronounced bubble coalescence higher up in the bed. Figure 8 shows a plot of the experimentally obtained bubble diameters at a constant excess gas velocity of 0.20 m/s and a jetsam mixture fraction of 75%. A total of over 500 points were plotted on the diagram. Predictions from Darton et al.’s37 semiempirical correlation are also plotted on the diagram. The black lines in the figure indicate the experimental values that are within (30% of the Darton et al.37 equation; the predictions from the semiempirical correlations are poorer higher up in the bed, where fewer experimental points are available. A general trend of an increase in the bubble diameter with an increase in the bed height can be discerned from the experimental scatter

Figure 9. Comparison of the simulated and experimental bubble diameter at an excess gas velocity of 0.20 m/s using a mixture fraction of 0.75 slag particles.

in Figure 8. Figure 9 shows a comparison between the simulated and experimental bubble diameters at different bed heights. The bubbles have been averaged into 10 classes by the method of interpolation. It can be observed that the experimental bubble diameter predictions are always larger than the simulated ones. A change in the mixture composition, for jetsam mixture fractions of 0.50 and 0.25, was found to have no effect on the numerical predictions of the bubble size in the computational simulations. A constant value for the bubble size was, in fact, observed for all other cases. This is because the same excess gas flow rate was applied in all of the simulations. Both simulations gave bubble size predictions that span a range of 0.005-0.085 m in concurrence with results reported earlier in Figure 7. Similarly, experimental observations showed that a change in the mixture composition of the jetsam had no noticeable effect on the experimental bubble size predictions. On the basis of the evidence reported in previous work (see Lettieri et al.39), a possible reason for the disparity between experimental and simulated diameters could be attributed to the dissimilarity in geometry between the 2D computational domain and the experimental “thin 3D” domain used in this work. Lettieri et al.39 demonstrated that the bubbles obtained from 2D simulations are always smaller that those obtained from 3D

6324 Ind. Eng. Chem. Res., Vol. 47, No. 16, 2008 Table 9. Experimental and Theoretical Values for Mixing Index M Obtained for Different Jetsam Fractions at an Excess Gas Velocity of 0.10 m/s mixing index M

coefficient of segregation

ωjet

CFD

exp

Wu and Baeyens29

CFD

exp

0.75 0.50 0.25

1.00 0.99 0.98

0.98 1.00 0.92

0.82 0.82 0.82

-0.304 +0.380 +0.736

-0.084 -0.065 +0.842

simulations, in agreement with experimental observation. The discrepancy between the CFD simulations and the results obtained from both experiments and the semiempirical correlation may also be attributed to the simplification in modeling of the distributor plate, as explained earlier. The omission of simulating the distributor plate may inherently introduce a deviation from predictions given by the semiempirical correlation, with the CFD bubble size prediction being smaller than the prediction by Darton et al. Furthermore, it is also worth noting that Darton et al. assumed in their paper that bubbles line up and grow as close together as possible, following specific paths, which is somewhat improbable in a chaotic system found in a fluidized bed. Hence, differences between experimental bubble size measurements and those obtained from the correlation may be due to the different bubble growth mechanisms that take place in experimental systems, i.e., the combination of the decrease in hydrostatic pressure toward the top of the fluidized bed, as well as bubbles coalescing horizontally and vertically through trailing bubbles and coalescence with neighboring bubbles, respectively. The above leads to the perpertuation of bubbles, experimentally, that might not necessarily line up and grow as close together as possible, as was assumed in the Darton model. 5.4. Segregation. For this study, the excess gas velocity, uxs, was reduced to 0.10 m/s after the computational bed was fully mixed, and simulations were carried out for a further period of 10 s in real time. This was done to enable the study of the segregation dynamics in the bed. The justification for not fluidizing at a superficial gas velocity slightly higher than the minimum fluidization velocity of the flotsam is based on the fact that segregation is promoted by the passage of bubbles in fluidized beds.36 Thus, neither mixing nor segregation will occur in a bubble-free bed at velocities slightly higher than the minimum fluidization velocity of the flotsam particle. Table 9 reports the mixing index M and the coefficient of segregation Cs calculated from experiments and after 10 s of simulation for all of the different case studies examined. The coefficient of segregation, introduced by Geldart,40 is given by Cs)

〈ωjet 〉b - 〈ωjet 〉t × 100 〈ωjet 〉b + 〈ωjet 〉t

(26)

where 〈ωjet〉b and 〈ωjet〉t are the mass fractions of the jetsam phase in the bottom and top halves of the bed. Clearly, Cs varies between -100 and +100, with -100 denoting perfect jetsam segregation at the top of the bed, 0 being representative of perfect mixing, and +100 indicating perfect jetsam segregation at the bottom of the bed. The values of M are also compared with the prediction obtained from the semiempirical correlation developed by Wu and Baeyens.29 All of the simulations show a good agreement with the experimental data and the semitheoretical expression albeit under circumstances of modest segregation. In qualitative terms, the numerical predictions, in concurrence with experimental results, show an increasing coefficient of segregation with a decrease in the jetsam fraction. This is due to the fact

that a decrease in the superficial gas velocity of the system, as a result of a decrease in the jetsam mass fraction, reduces the difference between the umf of the monocomponent jetsam particle (slag) and the superficial gas velocity applied. Quantitatively, with regards to the mixing index, the correlation by Wu and Baeyens29 yields predictions that deviate from experimental and computational measurements by 18%. The reason for this discrepancy might be that the mechanism of segregation is overwhelmed by the solids recirculation and exchange mixing mechanisms present in the experimental and computational beds. Copper and Coronella31 arrived at similar conclusions, following their numerical investigation into the segregation dynamics of natural rutile and coke particles using a commercial CFD package, FLUENT. In this case, they found that the mechanism of mixing pervaded in the fluidized bed over the antagonistic effect of segregation even when simulations were carried out at superficial gas velocities comparable to the minimum fluidization velocity of the jetsam particle. 6. Conclusions In this paper, an experimental validation of the numerical simulations of the mixing and segregation behavior of a bidisperse mixture of industrial rutile powders using a new model developed for bidisperse mixtures has been presented. Three different average compositions, corresponding to the average mass fractions of big particles (slag) of 0.25, 0.50, and 0.75 in the bed, were considered. The mixtures were fluidized at a constant excess gas velocity, uxs, of 0.20 m/s for 10 s. uxs was subsequently reduced to 0.10 m/s to enable a study of the segregation dynamics of the system. Qualitatively, the numerical simulations captured the mixing mechanism experimentally observed. A quantitative comparison between experimental and computational results on the mixing index expressed as a deviation percentage showed good agreement. Numerical results from macroscopic bed properties showed a good agreement with the experimental findings for bed properties like the bed height and voidage. However, a difference in the numerical prediction of the bubble holdup from experimental data was found. With regards of this, it was suggested that further investigations should be carried out on the effect of the solid stress tensor and, in particular, the role of the frictional viscosity of the mechanism of bubble growth in the fluidized bed. The effect of the distributor plate on the bubble formation and growth in the fluidized bed should also be investigated further. The different mechanisms that dominate bubble growth may be the underpinning factor to cause the scatter between the experimental, CFD, and predicted bubble size. In the segregation study, a modest segregation pattern was predicted both numerically and experimentally. The reason for this was attributed to the pervasiveness of the solids recirculation over the antagonistic effect of segregation in the CFD simulations, highlighting further the present gaps in current computer models in correctly predicting the mixing and segregation dynamics of these systems and the need for developing robust models capable of describing the particle-particle interaction effects in multiphase systems. While the materials in this paper have been simulated as a bidisperse mixture, work is currently ongoing to extend the model to polydisperse systems by solving the population balance model in conjunction with the conservation and momentum equations.

Ind. Eng. Chem. Res., Vol. 47, No. 16, 2008 6325

Acknowledgment The authors acknowledge financial support from the Engineering Physical Science Research Council and also from Tioxide Europe Ltd. In particular, many thanks are due to Roger Place and Dr. Stephen Sutcliffe. Nomenclature Roman Letters A ) area, m2 Ao ) catchment area, m2 Ar ) Archimedes number Cs ) coefficient of segregation D ) diameter, m db ) bubble diameter, m dR ) diameter ratio Ek ) elastic modulus F ) phase indicator G ) acceleration due to gravity, m2/s Ga ) Galileo number H ) height, m K ) phase indicator M ) mixing index mflot ) overall mass of flotsam, kg mjet ) overall mass of jetsam, kg P ) pressure, Pa QB ) volumetric flow rate, m3/s Re ) Reynolds number U ) superficial gas velocity, m/s umf ) minimum fluidization velocity, m/s uxs ) excess gas velocity, m/s Greek Letters Β ) fluid-solid interphase transfer coefficient, kg/m3 · s Ε ) fluid-phase volume fraction Σ ) solid-solid momentum transfer coefficient, kg/m3 · s Κ ) dilatational viscosity, Pa · s µf ) fluid shear viscosity, Pa · s Ff ) density of fluid, kg/m3 Fm ) density of the mixture, kg/m3 Φ ) solid volume fraction Ω ) particle mass fraction j ) stress tensor, N/m2 T Superscripts f ) fluid phase Fk ) phase Indictor K ) phase 1 or phase 2 p ) particle Subscripts 1, 2 ) phase indicator jet ) jetsam flot ) flotsam m ) mixture

Literature Cited (1) Owoyemi, O.; Mazzei, L.; Lettieri, P. CFD modelling of bi-disperse fluidized suspensions and investigation of the role of particle-particle drag on mixing and segregation. AIChE J. 2007, 53 (8), 1924. (2) Davidson, D. The role of computational fluid dynamics in process industries. Bridge 2002, 32, 4. (3) Gidaspow, D.; Ettehadieh, B. Fluidization in two dimensional beds with a jet π-hydrodynamic modelling. Ind. Eng. Chem. Fundam. 1983, 22, 193. (4) Gidaspow, D. Multiphase Flow and Fluidization; Academic Press Ltd.: London, 1994. (5) Jackson, R. Locally averaged equations of motion for a mixture of identical spherical particles and a Newtonian fluid. Chem. Eng. Sci. 1997, 52, 2457.

(6) Curtis, J.; van Wachem, B. Modeling particle-laden flows: A research outlook. AIChE J. 2004, 50, No 11. (7) Savage, S. Granular flows down rough inclinessreview and extension. In U.S.-Japan Seminar on New Models and ConstitutiVe Relations in the Mechanics of Granular Materials; Jenkins, J., Satake, M., Eds.; Elsevier: New York, 1982; p 261. (8) Johnson, P.; Jackson, R. Frictional-collisional consitutive relations for granular materials with application to plane shearing. J. Fluid Mech. 1987, 176, 67. (9) Foscolo, P. U.; Gibilaro, L. Fluid dynamic stability of fluidized suspensions: The Particle Bed Model. Chem. Eng. Sci. 1987, 42, 1489. (10) Makkawi, Y.; Ocone, R. A model for gas-solid flow in a horizontal duct with a smooth merge of rapid-intermediate-dense flows. Chem. Eng. Sci. 2006, 61, 4271. (11) van Wachem, B. G. M.; Schouten, J. C.; van den Bleek, C. M.; Krishna, R.; Sinclair, J. L. Comparative analysis of CFD models of dense gas-solid systems. AIChE J. 2001, 5, 1035. (12) Kuipers, J.; van Swaaji, W. Application of computational fluid dynamics to chemical reaction engineering. ReV. Chem. Eng. 1997, 13, 1. (13) van den Akker, H. On status and merits of computational fluid dynamics. In Proceedings of the 4th International Conference on Bioreactor and Bioprocess Fluid Dynamics; Nienow, W., Ed.; Edinburgh, Scotland, July 1997, BHR, ISBN 1-86058-090-4. (14) Sundaresan, S. Modelling the hydrodynamics of multiphase flow reactors: current status and challenges. AIChE J. 2000, 46, 1102. (15) Arastoopour, H. Numerical simulations and experimental analysis of gas/solid flow systems: 1999 fluor-daniel plenary lecture. Powder Technol. 2001, 119, 59. (16) Campbell, C. Granular material flowssan overview. Powder Technol. 2006, 162, 208. (17) Owoyemi, O.; Lettieri, P.; Place, R. Experimental validation of eulerian-eulerian simulations of rutile industrial powders. Ind. Eng. Chem. Res. 2005, 9996. (18) Owoyemi, O.; Lettieri, P. Stability analysis and CFD validation of a new fluid-particle interaction force for mono-component gas-solid fluidized beds. Powder Technol. 2008, 183 (1), 27. (19) Wallis G. B. One dimensional two-phase flow; McGraw-Hill: New York, 1969. (20) Gidaspow, D.; Syamlal, M.; Seo, Y. Hydrodynamics of fluidization of single and bi-disperse size particles: supercomputer modelling. Proc. 5th Eng. Found. Conf. Fluid. 1985, 5, 1. (21) Syamlal, M. The particle-particle drag term in a multiparticle model of fluidization; Topical Report; National Technical Information Service: Springfield, MA, 1987. (22) Bell, R. A. Numerical Modelling of Multi-Particle Flows in Bubbling Gas-Solid Fluidized Beds. Licentiate Thesis, Swinburne University of Technology, Melbourne, Australia, 2000. (23) Di Felice, R. Liquid Fluidization of Bi-disperse Solid Mixtures. Licentiate Thesis, University College London, London, 1988. (24) Gibilaro, L. Fluidization Dynamics; Butterworth Heinemann: London, 2001. (25) Mazzei, L.; Lettieri, P.; Elson, T.; Colman, D. A revised monodimensional particle bed model for fluidized beds. Chem. Eng. Sci. 2006, 61, 1958. (26) Gidaspow, D.; Syamlal, M.; Seo, Y. Hydrodynamics of fluidization of single and bi-disperse size particles: supercomputer modeling. Proc. 5th Eng. Found. Conf. Fluid. 1985, 5, 5. (27) Srinivasan, M. G.; Doss, E. D. Momentum transfer due to particleparticle interaction in dilute gas-solid flows. Chem. Eng. Sci. 1985, 40, 1791. (28) Gera, D.; Syamlal, M.; O’Brien, T. J. Hydrodynamics of particle segregation in fluidized beds. Int. J. Multiphase Flow 2004, 30, 419. (29) Wu, S. Y.; Baeyens, J. Segregation by size difference in gas fluidized beds. Powder Technol. 1998, 98, 139. (30) Davidson, J. F.; Harrison, D. Fluidised Particles; Cambridge University Press: Cambridge, U.K., 1963. (31) Gossens, W. R. A.; Dumont, G. L.; Spaepen, G. L. Fluidization of bi-disperse mixtures in the laminar flow region. Chem. Eng. Prog. Symp. Ser. 1971, 67, 38. (32) Cooper, S.; Coronella, C. J. CFD simulations of particle mixing in a bi-disperse fluidized bed. Powder Technol. 2005, 151, 27. (33) Anderson, J. Computational Fluid DynamicssThe basics with applications; McGraw-Hill: New York,1995. (34) Rowe, P. N.; Nienow, A. W.; Agbim, A. J. A preliminary quantitative study of particle segregation in gas fluidized bedssBidisperse systems of near spherical particles. Trans. Inst. Chen. Eng. 1972, 50, 324.

6326 Ind. Eng. Chem. Res., Vol. 47, No. 16, 2008 (35) Nienow, A. W.; Cheesman, D. J. The effect of shape on the mixing and segregation of large particles in gas-fluidized beds of small ones. Fluidization; Plenum: New York, 1980; Vol. 1, p 373. (36) Nienow, A.; Chiba, T. Fluidisation of dissimilar materials. Fluidization, 2nd ed.; Academic Press: London, 1985; Chapter 10. (37) Darton, R.; LaNauze, R.; Davidson, J.; Harrison, D. Bubble growth due to coalescence in fluidised beds. Trans. Inst. Chem. Eng. 1977, 55, 274. (38) Mazzei, L.; Lettieri, P. A numerical algorithm for the analysis of the bubble dynamics in 2D fluidized beds simulated by means of CFD multiphase codes. Int. J. Chem. React. Eng. 2006, 4, 1.

(39) Cammarata, L.; Lettieri, P.; Micale, G.; Colman, D. 2D and 3D CFD simulations of bubbling fluidized beds using Eulerian-Eulerian models. Int. J. Chem. React. Eng. 2003, 1, A48. (40) Geldart, D. Gas Fluidization Technology; John Wiley & Sons Ltd.: London, 1973.

ReceiVed for reView January 25, 2008 ReVised manuscript receiVed April 24, 2008 Accepted May 7, 2008 IE800138N