2390
Ind. Eng. Chem. Res. 2003, 42, 2390-2398
Hydrodynamic Simulations of Gas-Solid Flow in a Riser Lu Huilin*,† and Dimitri Gidaspow‡ Department of Power Engineering, Harbin Institute of Technology, Harbin 150001, China, and Department of Chemical and Environmental Engineering, Illinois Institute of Technology, Chicago, Illinois 60616
Conservation of mass and momentum equations for gas and solid phases were used to compute the hydrodynamics of flow in risers for two core-annular flow regimes. The principal input into the model was a semiempirical viscosity for fluidized catalytic cracking particles. The computed time-averaged radial fluxes and particle concentrations agree with experimental data obtained in a riser by Miller (Ph.D. Thesis, Illinois Institute of Technology, Chicago, IL, 1991) and Knowlton et al. (PSRI Challenge Problem Presented at the Eighth International Fluidization Conference, Tour, France, 1995). The model also computed the dynamics of flow, the dominant frequency of instantaneous particle concentrations, the granular temperature, and the granular pressure in reasonable agreement with measurements. 1. Introduction Circulating fluidized bed (CFB) hydrodynamic models are useful in understanding gas-solid mixing, scaleup, and control. The hydrodynamics of such a system impact the reactor performance: conversion, selectivity, and heat transfer. Furthermore, CFB riser operating conditions affect the efficiency of downstream equipment such as cyclones, filters, standpipes, and so forth. Several modeling efforts employing very different mathematical formulations have appeared in the literature to predict the relationship between the solid concentration, operating conditions, and riser geometry. The models can be classified into three broad categories:1,2 (I) those that predict the axial solid suspension density profile but not the radial profile; (II) those that predict the radial profile by assuming two or more regions, such as core-annulus or clustering annular flow models; (III) those that employ the fundamental equations of fluid dynamics to predict two-phase gas-solid flow. Type III models, because of their generality, are suitable for investigating the riser local flow structure and predicting the effects of complex geometry, such as corner effects in CFB combustors. Type III models are based on the principles of conservation of mass, momentum, and energy for each phase. The early models were developed by Davidson,3 Jackson,4 and Soo.5 Gidaspow6 has critically reviewed the models. Two types of hydrodynamic models were used to predict riser flow: the viscosity model and the kinetic theory model with and without gas turbulence. The kinetic theory model was applied to model riser flow.7-12 The kinetic theory approach uses one equation model to determine the turbulent kinetic energy of the particles (granular temperature) and allows the determination of the pressure and viscosity of the solid in place of empirical relations. Tsuo and Gidaspow13 were among the first to predict the core-annular flow regime using the viscosity model. They successfully used an empirical relation for the solid viscosity that was later used by Sun and Gidaspow.14 Type III models form the * To whom correspondence should be addressed. E-mail:
[email protected]. † Harbin Institute of Technology. ‡ Illinois Institute of Technology.
theoretical basis of the recently formed Multiphase Fluid Dynamics Research Consortium.15 It consists of American chemical companies, led by Dow Chemical, six U.S. national laboratories funded by the U.S. Department of Energy, Office of Industrial Technology, and six universities, supported by industry. Its mission is to solve the gas-solid riser problem for Geldart type A and C particles. As part of the work of these efforts, we show in this paper that the viscosity model is able to describe the most significant feature of gas-solid twophase flow hydrodynamics and the time-averaged coreannular flow regime in the CFBs. 2. Model Input Table 1 summarizes the equations that are well-posed as an initial value problem and therefore can be solved by stable finite difference techniques. Arastoopour and Gidaspow16 established four different two-phase flow models for the description of a one-dimensional steadystate pneumatic conveying system. Tsuo and Gidaspow13 were the first to simulate cluster formation and predicted the core-annular flow of particles by using a two-dimensional transient two-phase flow model with a constant viscosity of particles as an input. Arastoopour et al.8 similarly developed a steady-state two-dimensional gas-solid flow model using the method of lines. The model with solids viscosity as an input was also used by Benyahia et al.17 to model gas-solid flow in a riser. Sun and Gidaspow14 predicted gas-solid behavior to input an empirical equation of particle viscosity proposed by Miller and Gidaspow18 into the model. The viscosity of particles was used to model the viscous effects of the solid phase. The viscosity input model was also used reasonably successfully to model flow patterns of bubbles in bubbling fluidized beds.19-24 The model presented in Table 1 has two key empirical inputs which come from the constitutive equation for the particles: the particle pressure and viscosity. The particle pressure, Ps, is mainly caused by particle interaction (collisions) in the dense gas-solid mixture flow. This pressure transmits a force both by the shortduration collisional impacts and by the long-duration particle-particle contact. The effect of the particle pressure-gradient term is to keep the particles apart so that the calculated particle concentration does not
10.1021/ie020521q CCC: $25.00 © 2003 American Chemical Society Published on Web 04/25/2003
Ind. Eng. Chem. Res., Vol. 42, No. 11, 2003 2391 Table 1. Mathematical Model of Gas-Solid Flow A. conservation laws (1) continuity equations for different phases (a) fluid phase
∂ ( F ) + ∇(gFgvg) ) 0 ∂t g g
(3)
(b) particulate phase
∂ ( F ) + ∇(sFsvs) ) 0 ∂t s s
(4)
(2) momentum equations (a) fluid phase
∂ ( F v ) + ∇(gFgvgvg) ) ∂t g g g -∇Pg + gFgg + β(vs - vg) + ∇[τg] (5) (b) particulate phase
∂ ( F v ) + ∇(sFsvsvs) ) -∇Ps + sFsg + β(vg - vs) + ∇[τs] (6) ∂t s s s Figure 1. Variations of the solid viscosity.
exceed the maximum concentration obtainable for a given size distribution of the particles. Several models for the particle pressure-gradient term presented in the literature are based on the following equation: ∇Ps ) -G(g) ∇g. The function G(g) can be thought of as a modulus of elasticity for the particulate phase, which must be modeled empirically. An overview of different models for the modulus of elasticity is given by Massoudi et al.25 Equation 10 of the modulus of elasticity was used6 in this study. The results obtained from the different models of G(g) differ widely. However, the main effect of the particle pressure-gradient term is only to prevent the particulate phase from becoming too densely packed, and the deviations between the different modes for G(g) may not significantly affect the timeaveraged solutions for the concentration, velocity, and pressure fields of the phases. This remains to be investigated in the future. The viscous stress of the particle phase enters into the momentum equation. The stress tensor of particles can be taken into account by modeling the particle viscosity. For the measured particle concentration by means of a X-ray densitometer, solid mass flux measured using an extraction probe, and calculated particle velocity distribution from mass conservation of particles of 75 µm FCC particles in a riser, Miller and Gidaspow18 determined the solid viscosity from a mixture momentum balance, neglecting transient effects and assuming that the gas and solid velocity gradients are of the same order. They found that the solid viscosity can be approximated by the solid concentration multiplied by a constant: µs ) 5.0s (P). More recent particle viscosity data of FCC particles obtained by Gidaspow et al.26,27 using the particle image velocity technique (see Figure 1) show that a better correlation for viscosity, partially based on the kinetic theory of granular flow, is expressed as follows:
µs ) 0.165s1/3g0
(1)
where g0 is the radial distribution function at contact. Kinetic theory of granular flow shows that the granular temperature rises as the two-third power of the solid volume fraction for dilute flow. This is analogous to the rise of thermal temperature upon compression, where 2/3 is the ratio of specific heats. Because
kinematic viscosity ) mean free path × square root of the granular temperature (2) where the square root of the granular temperature is
B. constitutive equations for phases (1) constitutive equations for stresses (a) fluid-phase stress
2 [τg] ) µgg[∇vg + ∇vgT] - µgg∇vg 3
(7)
(b) particulate-phase stress
2 τs ) µss[∇vs + ∇vsT] - µss∇vs 3
(8)
(2) empirical particulate-phase pressure and viscosity models (a) particulate pressure
∇Ps ) G(g) ∇s G(g) ) 10
-8.686g+6.385
(9) (10)
(b) FCC particle viscosity
µs ) 0.165s1/3g0 (P) 1/3 -1
g0 ) [1 - (s/s,max) ]
(11) (12)
(3) fluid-phase viscosity model (SGS model)
µg ) µg,l + Fg(0.1∆)2(τg‚τg) with ∆ ) x∆r∆z (4) fluid-particulate interphase drag coefficients 3 Fgk|vg - vk| -2.65 g g 0.8: βg,s ) CD g 4 d 2 150s µg Fgs|vg - vs| g < 0.8: βg,s ) + 1.75 2 gd (gd) 24 CD ) (1 + 0.15Re0.687) Re < 1000 Re Re g 1000 CD ) 0.44 Fgg|vg - vs|d Re ) µg (5) ideal gas equation Fg ) p/RT
(13)
(14) (15) (16) (17) (18) (19)
the random particle oscillation velocity, the one-third dependence of the volume fraction is obtained for dilute flow by elimination of the granular temperature. The radial distribution function accounts for the probability of particle collisions and is near 1 when the flow is dilute and becomes infinity when the flow is so dense that motion of particles is impossible. It can be expressed in terms of the equation by Carnahan and Starling28 or Bagnold.29 Experimental data for FCC particles was found to lie between these theoretical expressions.30 At high volume fractions, both expressions give large values, raising the viscosity rapidly near maximum packing. This paper presents a model of an isothermal gassolid two-phase flow. In the solid phase, a semiempirical particulate viscosity is incorporated. Radial solids flux and particle concentration distributions are computed using this improved viscosity correlation. The model is verified against experimental data of the CFB with
2392
Ind. Eng. Chem. Res., Vol. 42, No. 11, 2003
different gas velocities and solid mass fluxes by Miller31 and Knowlton et al.32 3. Mathematical Model and Solution Method The equations used for the prediction are the conservation of mass and momentum for the solid and gas, given in Gidaspow’s book6 as model B. They are summarized in Table 1. With a substitution of a Newtonian-type stress for each phase, the differential equations (3)-(6) express conservation principles. The differential equations express the conservation over an infinitesimal control volume and need to be discretized over finite volumes.33 The simulations were carried out with a modified CFD code K-FIX, which was previously used to model the flow in the fluidized beds (ref 6, pp 216-228). This software allows free implementation of extra equations, boundary conditions, and differencing schemes. For solving the difference equations obtained from the differential equations, the higher order total variation diminishing (TVD) scheme is used. This TVD scheme incorporates a modification to the higher order upwind scheme. The solution of the pressure from the momentum equations requires a pressure correction equation, correcting the pressure and the velocities after each iteration of the discretized momentum equations. The calculated pressure is used to compute the density of the gas phase. For the gas phase, the no-slip boundary condition at the walls is used:
ug,w ) 0, vg,w ) 0
Figure 2. Schematic drawing of a two-dimensional riser with inlet and initial conditions.
(20)
For the solid phase, particles can be moved while touching the wall. In this work the free slip boundary condition is chosen.34 The vertical and radial velocities of particles at the wall are calculated by
us,w ) 0, vs,w ) -
d ∂vs 1/3 ∂n
(21)
s
where the n direction is normal to the wall. At the bottom of the riser, the gas inflow is specified. The pressures in the mesh cells at the top of the riser are fixed at a specific value. In this study, the simulations are carried out over a two-dimensional plane in a cylindrical coordinate system. The common problem encountered by using such coordinates is the singularity that occurs in solving above differential equations when the radial distance r equals zero. This problem can be overcome by carefully setting up the simulation meshes. The number of cells was chosen so that the radial distance r ) 0.0 never occurs as a denominator in the finite-differenced equations. 4. Computational Results Knowlton et al.32 presented the experimental data obtained from a CFB with FCC particles. The unit consists of a 20 cm i.d. riser with a length of 14.5 m. The measurement was located at the height of 3.9 m from the inlet. Figure 2 shows the riser section used in the present numerical simulation of gas-solid flow. The computational domain is considered to be of 23 grids radially and 150 nonuniform grids along the axis of the riser. A total of 3450 fluid cells resulted from this grid distribution. A constant time step of 5 × 10-5 s was
Figure 3. Solid velocity vector plots at 20, 25, and 29 s of real simulation time.
used. Figure 3 shows the simulated flow pattern of the transient particle velocity at the gas velocity and mass flux of 7.6 m/s and 489 kg/m2‚s, respectively. Distributions refer to three different times (20, 25, and 29 s) from the beginning of the simulation. A characteristic feature of the flow is the oscillating motion of solid particles from one wall to the other wall through the center line of the riser. The complex and transient velocity fields of particles flowing up and down are evident. It should be noted how downflow conditions at one wall are mostly associated with an upflow condition at the other wall, but more complex combinations are possible too. To compare simulation results with Knowlton et al. data,32 time-averaged distributions of flow variables have been computed. Parts a and b of Figure 4 show the computed time-averaged solid mass flux distributions at gas velocities of 5.2 and 7.6 m/s for the same solid mass flux of 489 kg/m2‚s. In these cases, to obtain the correct time average, the simulation time continued
Ind. Eng. Chem. Res., Vol. 42, No. 11, 2003 2393
Figure 4. Simulated solid mass flux distributions.
Figure 5. Simulated solid mass flux distributions.
to 60 s of real time. From Figure 4a, it can be seen that the parabolic distribution of mass flux is predicted with downflow at the wall for the case of the gas velocity and solid mass flux of 5.2 m/s and 489 kg/m2‚s. Although the values have quantitative discrepancies between them, the trends are the same. As the gas velocity increases to 7.6 m/s at the same solid mass flux, the downflow of particles near the wall disappeared. The mass solid flux decreases toward to the wall. However, at the center of the riser, the simulation results obtained using the equation of Miller and Gidaspow and eq 1 reach a local minimum point, similarly to experimental data. From these figures, we see that the simulation results predicted using the equation of Miller and Gidaspow18 are different from those using eq 1. At the center region of the riser, the values obtained using the equation of Miller and Gidaspow are smaller than those using eq 1 and experimental data. The simulated results predicted using the equation of Miller and Gidaspow, however, are larger than those using eq 1 and experiments at the wall region. We see that the computed results obtained using eq 1 agree well with experimental data. Parts a and b of Figure 5 show the simulated timeaveraged solid mass flux distributions at solid mass fluxes of 489 and 189 kg/m2‚s for the same gas velocity of 11 m/s. From these figures, we see that the upflow of particles occurs at high gas velocity and relatively low solid mass flux. The maxima in both the simulations and the experimental data are between the wall and the center of the riser. At the center region of the riser,
the values obtained using the equation of Miller and Gidaspow18 are larger than those using eq 1. The simulated results predicted using the equation of Miller and Gidaspow, however, are smaller than those using eq 1 at the wall region. The simulations obtained using eq 1 captured the uncommon mass flux behavior at the higher gas velocities and solid mass fluxes in the riser. Parts a and b of Figure 6 show the power spectrum density profiles based on the computed instantaneous solid concentrations at the location of 8.5 m from the riser inlet by means of the fast Fourier transform method for the gas velocity and solid mass flux of 7.6 m/s and 489 kg/m2‚s and 11 m/s and 189 kg/m2‚s, respectively. The results predicted using the equation of Miller and Gidaspow18 are also indicated. The main reason for calculating the frequency of oscillations of gas-solid flow was to know the minimum time required to conduct proper time averaging to evaluate the granular temperature of particles discussed below. These plots give the dominant frequency that is the highest values existing in the frequency domain. From these figures we see that the dominant frequency predicted using the equation of Miller and Gidaspow was lower than that using eq 1. The dominant frequencies obtained using eq 1 are from about 0.07 Hz at the gas velocity and solid mass flux of 7.6 m/s and 489 kg/ m2‚s to 0.01 Hz at the gas velocity and solid mass flux of 11 m/s and 189 kg/m2‚s. These values are compared with IIT CFB data based on the measurements of the solid concentration by means of a γ-ray densitometer.35 Figure 7 shows a comparison of the computed dominant
2394
Ind. Eng. Chem. Res., Vol. 42, No. 11, 2003
Figure 6. Profile of the power spectrum density.
Figure 7. Comparison of the dominant frequency to the analytical solution.
Figure 8. Granular temperature as a function of the solid volume fraction in any cell of the computational domain (the inset shows the experimental granular temperature values from ref 30).
frequency to the experimental data. The theory of fluidization shows that the mechanism of large-scale oscillations is a coupled pressure and density nonlinear wave propagation phenomenon (ref 6, pp 131-134). For a dense gas-solid flow, the density and the flux obey a simple wave equation with densities propagating at pseudosonic velocities of particles, Cs. For a dilute flow, Gidaspow et al.36 give the expression of the frequency of oscillation of particles, f (Hz), as
f)
()[
1 g 2π x0
1/2
]
(1 + 4.65s/g)s s,0
1/2
(22)
where s,0 and x0 are some initial solid volume fraction and height filled with particles. We see that the basic frequency is that caused by gravity, (g/x0)1/2, and that the frequency becomes very small as the volume fraction of the particles becomes small. The frequencies predicted by eq 22 are given in Figure 7. Very low frequencies were also obtained by Sun and Gidaspow.14 Equation 22 suggested that for very dilute flow the frequency will vanish, giving a steady state. In such a situation, the frequency will jump to a value of the gas flow alone. From the computed instantaneous velocities of particles, the granular temperature is calculated, assuming the angular and radial velocity distributions to be equal. Experimental measurements by means of the particle image velocity technique show that the vertical component of the particle velocity is larger than that of the
Figure 9. Solid pressure as a function of the solid volume fraction in any cell of the computational domain (the inset shows the experimental solid pressure values from ref 30).
radial component in a CFB riser.27 Although the particle velocities of the angular and radial components are not the same, to a first approximation, they are assumed to be equal. The granular temperature of the particles, θ, is calculated as follows:
1 2 θ ) σz2 + σr2 3 3
(23)
where σz and σr are the axial and radial standard
Ind. Eng. Chem. Res., Vol. 42, No. 11, 2003 2395
Figure 10. Distribution of the particle density.
deviations of the particle velocity. The standard deviation σ is given by
σ)
x
1
N
(vi - Vs)2 ∑ N - 1i)1
(24)
where vi is the instantaneous particle velocity and Vs the mean value of the particle velocity. Figure 8 shows the behavior of the granular temperature, θ, as a function of the solid volume fraction in any cell of the riser. The box included illustrates the same dependence as results from the experimental data.30 Similar plots can be computed at different conditions but are not reported here for the sake of brevity. It can be seen that the relation between the granular temperature and the solid volume fraction is not universal. However, the diagram presents some interesting features that confirm the experimental evidence. First, the granular temperature exhibits a maximum near a solid volume fraction at about 0.05-0.10. A similar maximum is expected at a solid volume fraction of about 0.05-0.10 based on experiments near the wall. Second, the granular temperature tends to decrease with the solid volume fraction at the high concentration of solids. Finally, at the dilute conditions the granular temperature increases with an increase of the solid volume fraction. Such a dependence in the dilute conditions was proven by Gidaspow and Huilin.30 A similar behavior of the granular temperature vs solid volume fraction is also described by Benyahia et al.17 It can be seen that the predicated granular temperature is smaller than that of experiments. The most probable reason could be that the large-scale oscillations that were measured experimentally were not included in the computational results. The time-averaged experimental measurements always include small- and large- scale oscillations, and the computational results include only the small-scale oscillations. Therefore, a direct comparison between the computation and experiment would be accurate only if the large-scale oscillations are subtracted from the experimental data. Based on the kinetic theory of granular flow, the solid pressure, Ps, can be predicted as a function of the solid volume fraction and granular temperature as follows:6
Ps ) sFsθ[1 + 2(1 + e)g0s]
(25)
Figure 9 shows the computed solid pressure as a function of the solid volume fraction in any cell of the
Figure 11. Distribution of the solid concentration at 32 and 39 s.
riser. The box included illustrates the results from the experimental data by Gidaspow and Huilin.30 The diagram clearly shows the difference between the computation and experiment. This was mainly due to the underprediction of the granular temperature. However, the trend is the same. The solid pressure tends to increase with an increase of the solid volume fraction. Parts a and b of Figure 10 show the computed and experimental solid density profiles at the gas velocity and solid mass flux of 5.2 m/s and 489 kg/m2‚s and 11 m/s and 196 kg/m2‚s, respectively. All of the figures show a core-annular flow behavior. In the center of the riser, the concentration of particles is low. It is high at the wall. The particle concentration segregation is mainly due to the stagnant wall region and motion of the particles to the wall region by oscillations. There is a reasonable agreement between the computation and experiment. Miller31 has studied the hydrodynamics of gas-solid flow in the riser. He obtained the radial profile data at different heights in the riser. The X-ray densitometer was used to measure the particle concentration, and the solid flux was measured by means of an extraction probe. The riser is 6.6 m high with an i.d. of 0.075 m. The mean diameter and density of FCC particles are
2396
Ind. Eng. Chem. Res., Vol. 42, No. 11, 2003
Figure 12. Experimental and simulated particle concentration profiles.
Figure 13. Particle velocities as a function of dimensionless distance.
Figure 14. Effect the boundary conditions of the particle phase on simulated particle velocities.
75 µm and 1650 kg/m2‚s, respectively. Figure 11 shows the simulated flow pattern of the transient particle concentration at the gas velocity and mass flux of 2.89 m/s and 20.4 kg/m2‚s, respectively. Distributions refer to two different times (32 and 39 s) from the beginning of the simulation. The behavior of the flow is characterized by the formation of particle clusters near the wall. Such portions of particles undergo a vigorous up and down motion, thus favoring a particle recirculation all over the riser. The nonhomogeneities of the particle concentration are evident. Parts a and b of Figure 12 show the profile of the time-averaged particle concentration at the gas velocity and solid mass flux of 2.89 m/s and 12.0 kg/m2‚s, respectively. The simulated
results obtained using the equation of Miller and Gidaspow18 are also shown in the figures. The computed particle concentrations obtained using the equation of Miller and Gidaspow and eq 1 have the same trends. Comparing the simulation with the experimental data, we can see that the particle concentration obtained using the equation of Miller and Gidaspopw near the walls is lower than that of experimental data. The coreannular flow structure in the riser is observed. From these figures, it is found that there is a reasonable agreement between the simulation results obtained using eq 1 and those of the experiments. Parts a and b of Figure 13 show the profiles of the particle velocity at the gas velocity and solid mass flux
Ind. Eng. Chem. Res., Vol. 42, No. 11, 2003 2397
of 2.89 m/s and 20.4 kg/m2‚s, respectively. At the height of 5.52 m, the predicted particle velocity obtained using the equation of Miller and Gidaspow is lower than experimental data. The simulation results predicted using eq 1, however, are in agreement with experimental data. At the height of 4.46 m, we see that the particle velocities obtained using the equation of Miller and Gidaspow18 and eq 1 are higher than experimental data in the center. The values predicted using the equation of Miller and Gidaspow are underpredicted at the walls. From these figures, it can be seen that the particle velocity is negative near the wall, and it is positive at the center of the pipe. This means that the particles flow up at the center by gas forces and flow down near the wall by the gravity force of the particles. It can be seen that there is an agreement between the Miller data and our simulations. Boundary conditions for the particulate phase at the wall are much more complicated than those for the gas phase. One may assume that the no-slip boundary condition is imposed for the particulate phase at the wall; that is, the parallel and perpendicular velocity components of the particles are to be zero at the wall. Hence, the slip boundary condition of the particle phase at the wall is expressed as us,w ) 0.0 and vs,w ) 0.0. Parts a and b of Figure 14 show the distributions of the particle velocity at the gas velocity and solid mass flux of 3.48 m/s and 20.4 kg/m2‚s, respectively. For the assumption of the slip boundary condition, the predicted particle velocity decreases toward to the wall, reaches a minimum value near the wall, and then goes to zero at the wall. However, the predicted particle velocity obtained using the no-slip boundary condition decreases to zero at the walls. No backflow of the particles is experimentally measured. This trend is not captured using the assumption of the slip boundary condition of the particle phase at the walls. 5. Conclusions The observed core-annular flow structure in a riser was predicted by a model that solves the transient conservation of mass and momentum equations for the particle and gas phases. The input into the model was the measured particulate viscosity of FCC particles and an estimate of the solid modulus of elasticity. Frequencies of instantaneous particle concentrations are compared to an equation derived from the equations of motion. Using a semitheoretical equation of the solid viscosity of FCC particles, the predictions are close to experimental data of Knowlton et al.32 and Miller.31 For high solids flux, the model predicts two types of core-annular flow regimes: a regime with a parabolic flux and downflow at the wall and a regime with a low flux at the pipe center and a maximum near the wall with no downflow. The latter flow regime may reduce undesirable backmixing in reactors. The two-fluid model equations must be solved together with appropriate constitutive equations of gas and particle phases. The particle-phase turbulence can be modeled with the kinetic theory of granular flow. Alternatively, The simply way to implement turbulence in the particulate phase is to include the concentrationdependent particle viscosity, using a Newtonian stressstrain relation. Although the linear relations between the particle viscosity and particle concentration, such as the Miller and Gidaspow equation, are easy to implement into the code, its validity used in the high
concentration of gas-particle fluidization is, at least, questionable. The proposed eq 1 of particle viscosity depends on the radial distribution function that accounts for the probability of particle collisions. It can be anticipated that this equation can be extended to dense gas-solid two-phase flow. Acknowledgment This study was supported by the National Science Foundation in China through Grant 10072019 and the U.S.A. National Science Foundation through Grant CTS-0086250. Nomenclature c ) constant Cs ) sonic velocity of the particle d ) particle diameter D ) diameter of the riser e ) restitution coefficient Gs ) solid mass flux g ) gravity g0 ) radial distribution function H ) height I ) unit tensor Pg ) fluid pressure Ps ) particle pressure r ) radii from the centerline of the riser ∆r ) r direction width of the control volume R ) universal gas constant, radii of the riser T ) absolute temperature of the gas phase U ) velocity ug ) gas velocity V ) mean velocity vg ) gas velocity vs ) particle velocity x0 ) initial height filled with particles ∆z ) z direction width of the control volume Greek Letters τg ) gas stress tensor τs ) stress tensor of the particulate phase ξ ) bulk viscosity µs ) particulate viscosity g ) voidage s ) solid volume fraction s,max ) solid volume fraction at packing θ ) granular temperature σ ) standard deviation Fs ) particle density Fg ) gas density β ) gas-particle interphase drag coefficient Subscripts g ) gas phase s ) solid phase w ) wall
Literature Cited (1) Harris, B. J.; Davidson, J. F. Axial and radial variation of flow in circulating fluidized bed riser. In Circulating Fluidized Bed Technology; Avidan, A. A., Ed.; AIChE: New York, 1994; Vol. IV, pp 103-109. (2) Berruti, F.; Chaouki, J.; Godfroy, L.; Pugsley, T. S.; Patience, G. S. Hydrodynamics of circulating fluidized bed risers: a review. Can. J. Chem. Eng. 1995, 73, 579. (3) Davidson, J. F. Symposium on fluidization: Discussion. Trans. Inst. Chem. Eng. 1961, 39, 230.
2398
Ind. Eng. Chem. Res., Vol. 42, No. 11, 2003
(4) Jackson, R. The mechanics of fluidized bed: I the stability of the state of uniform fluidization. Trans. Inst. Chem. Eng. 1963, 41, 13. (5) Soo, S. L. Fluid dynamics of multiphase systems; BlaisdellGinn: Waltham, MA, 1967. (6) Gidaspow, D. Multiphase flow and fluidization: Continuum and kinetic theory descriptions; Academic Press Inc.: Boston, 1994. (7) Sinclair, J. L.; Jackson, R. Gas-particle flow in a vertical pipe with particle-particle interaction. AIChE J. 1989, 35, 1473. (8) Arastoopour, H.; Pakdel, P.; Adewumi, M. Hydrodynamic analysis of dilute gas-solid flow in a vertical pipe. Powder Technol. 1990, 62, 163. (9) Louge M.; Mastorakos, E.; Jenkins, J. The role of particle collisions in pneumatic transport. J. Fluid Mech. 1991, 231, 345. (10) Pita, J. A.; Sundareasan, S. Gas-solid flow in vertical tubes. AIChE J. 1991, 37, 1009. (11) Samuelsberg, A.; Hjertager, B. J. H. Computational modeling of gas/particle flow in a riser. AIChE J. 1996, 42, 1536. (12) Goldschmidt, M. J. V.; Kuipers, J. A. M.; van Swaaij, W. P. M. Hydrodynamic modelling of dense gas-fluidized beds using kinetic theory of granular flow. Chem. Eng. Sci. 2001, 56, 571. (13) Tsuo, Y. P.; Gidaspow, D. Computation of flow patterns in circulating fluidized beds. AIChE J. 1990, 36, 885. (14) Sun, B.; Gidaspow, D. Computation of circulating fluidized bed riser flow for the Fluidization VIII benchmark test. Ind. Eng. Chem. Res. 1999, 38, 787. (15) Thompson, T. B. Organizer, Dow Chemical. The Multiphase Fluid Dynamics Research Consortium; U.S. Department of Energy, Office of Industrial Technology: 1999. (16) Arastoopour, H.; Gidaspow, D. Vertical pneumatic conveying using four hydrodynamics Models. Ind. Eng. Chem. Fundam. 1979, 18, 123. (17) Benyahia, S.; Arastoopour, H.; Knowlton, T. M. Prediction of solids and gas flow behavior in a riser using a computational multiphase flow approach. In Fluidization; Fan, L. S., Knowlton, T. M., Eds.; Engineering Foundation: New York, 1998; Vol. IX, pp 493-500. (18) Miller, A.; Gidaspow, D. Dense, vertical gas-solid flow in a pipe. AIChE J. 1992, 38, 1801. (19) Kuipers, J. A. M.; van Duin, K. J.; van Beckun, F. P. H.; van Swaaij, W. P. M. A numerical model of gas fluidized beds. Chem. Eng. Sci. 1992, 47, 1913. (20) Lyczkowsky, R. W.; Gamwo, I. K.; Dobran, F.; Ali, H.; Chao, B. T.; Chao, M.; Gidaspow, D. Validation of computed solids hydrodynamics and pressure oscillation in bubbling atmospheric fluidized bed. Powder Technol. 1993, 76, 65. (21) Anderson, K.; Sundaresan, S.; Jackson, R. Instabilities and the formation of bubbles in fluidized beds. J. Fluid Mech. 1995, 303, 327.
(22) Enwald, H.; Peiarno, E.; Almstedt, A. E. Eulerian twophase flow theory applied to fluidization. Int. J. Multiphase Flow 1996, 22, 21. (23) Go¨z, M. F.; Sundaresan, S. Growth, saturation and scaling behavior of one- and two-dimensional disturbances in fluidized beds. J. Fluid Mech. 1998, 362, 83. (24) Gustavsson, M.; Almstedt, A. E. Numerical simulation of fluid dynamics in fluidized beds with horizontal heat exchange tubes. Chem. Eng. Sci. 1999, 55, 857. (25) Massoudi, M.; Rajagopal, K. R.; Ekmann, J. M.; Mathur, M. P. Remarks on the modeling of fluidized systems. AIChE J. 1992, 38, 471. (26) Gidaspow, D.; Huilin, L.; Mostofi, R.; Wu, Y. Turbulence, viscosity and numerical simulation of FCC particles in CFB. Fluidization and Fluid-particle Systems, AIChE Annual Meeting, Los Angeles, 1997; pp 58-62. (27) Gidaspow, D.; Huilin, L. Collisional viscosity of FCC particles in a CFB. AIChE J. 1996, 42, 2503. (28) Carnahan, N. F.; Starling, K. E. Equation of state for nonattracting rigid spheres. J. Chem. Phys. 1969, 51, 635. (29) Bagnold, R. A. Experiments on a gravity-free dispersion of large solid spheres in a Newtonian fluid under shear. Proc. R. Soc. London 1954, A225, 49. (30) Gidaspow, D.; Huilin, L. Equation of state and radial distribution function of FCC particles in a CFB. AIChE J. 1998, 44, 279. (31) Miller, A. Dense, vertical gas-solid flow in a pipe. Ph.D. Thesis, Illinois Institute of Technology, Chicago, IL, 1991. (32) Knowlton, T.; Geldart, D.; Masten, J.; King, D. Comparison of CFB hydrodynamic models. PSRI Challenge Problem Presented at the Eighth International Fluidization Conference, Tour, France, May 1995. (33) Patankar, S. V. Numerical Heat Transfer and Fluid Flow; Hemisphere: New York, 1980. (34) Ding, J.; Gidaspow, D. A bubbling fluidization model using kinetic theory of granular flow. AIChE J. 1990, 36, 523. (35) Gidaspow, D.; Huilin, L.; Therdthianwong, A. Measurement and computation of turbulence in a circulating fluidized bed. Fluidization, Tour, France (preprints), 1995; Vol. VIII; pp 81-88. (36) Gidaspow, D.; Huilin, L.; Mostofi, R. Large scale oscillations or gravity waves in riser and bubbling beds. In Fluidization; Kwauk, M., Li, J., Yang, W. C., Eds.; Engineering Foundation: New York, 2001; Vol. X, pp 317-324.
Received for review July 15, 2002 Revised manuscript received January 22, 2003 Accepted January 30, 2003 IE020521Q
