Novel Approach to Characterize Fluidized Bed Dynamics Combining

These experimental data are exported to a finite element software where the pressure and gas velocity fields are obtained numerically. ... The effect ...
0 downloads 3 Views 14MB Size
5010

Ind. Eng. Chem. Res. 2009, 48, 5010–5023

Novel Approach to Characterize Fluidized Bed Dynamics Combining Particle Image Velocimetry and Finite Element Method J. A. Almendros-Iba´n˜ez,*,† D. Pallare`s,‡ F. Johnsson,‡ and D. Santana§ Escuela de Ingenieros Industriales, Dpto. de Meca´nica aplicada e Ingenierı´a de Proyectos, Castilla La Mancha UniVersity, Campus UniVersitario s/n, 02071, Albacete, Spain, Renewable Energy Research Institute, AVda. de la InVestigacio´n s/n, 02071, Albacete, Spain, Chalmers UniVersity of Technology, Department of Energy and EnVironment, SE 412 96 Go¨teborg, Sweden, and UniVersidad Carlos III de Madrid, ISE Research Group, Thermal and Fluid Engineering Department, AVda. de la UniVersidad 30, 28911, Legane´s, Madrid, Spain

This work presents a new experimental-numerical method combining PIV (particle image velocimetry), digital image analysis, and FEM (finite element method) to obtain gas and particle motion around bubbles in a 2D freely bubbling fluidized bed. The bubble geometry is captured with a high speed video camera while the particle velocity is measured using a PIV technique. These experimental data are exported to a finite element software where the pressure and gas velocity fields are obtained numerically. The flow equations proposed by Davidson’s model have been chosen to exemplify the application of the method presented in this paper. Different bubble types have been analyzed: slow and fast bubbles but also erupting and interacting bubbles. The effect on the gas flow of bubbles with a nonzero horizontal velocity component has also been analyzed, and it is shown how such bubbles interchange gas with the main stream. In addition, the Darcy’s law included in Davidson’s model has been extended including a quadratic term (i.e., Ergun’s equation) to take into account non-Darcy effects and the PIV results have been used to evaluate the error committed when voidage variation around bubbles is neglected. The results obtained including non-Darcy effects show differences in the gas velocity magnitude in regions near the nose and the wake of the bubble, where the gas velocity and, consequently, the local Reynolds number are higher. Nevertheless, these local differences give no significante difference in the gas streamlines around bubbles. Finally, the study of the voidage variation around bubbles shows that this variation is important only in a region of small thickness adjacent to the bubble. Introduction Fluidized bed technology has various applications in different fields of the industry. Some examples are the following: fluid catalytic cracking (FCC), biomass and coal gasifiers, and combustors or dryers. In most of the industrial applications (except FCC) the fluidizing powders used as inert particles to maintain a good fluidization quality are type B according to Geldart’s classification.1 For this type of particle, the minimum fluidization velocity is approximately equal to the minimum bubbling velocity. Therefore, when the bed is fluidized, bubbles appear, grow, and coalesce along the bed height until they reach the bed surface, i.e. these bubbles determine the pressure distribution as well as the dominant frequency in the bed.2,3 The bubbles have also an influence on the particle dispersion and therefore on the mixing in the bed.4 In the beginning of the 1960s, different theoretical models were developed to model particle and gas flow around bubbles.5-9 From the 1980s to the present, there has been a growing interest in the use of computational fluid dynamics (CFD) simulations in order to model the behavior of a fluidized bed. The development of the CFD simulation tool also benefits from the continuous increase in computer performance. van Wachem and Almstedt10 reviewed the state of the art in this field. Two main approaches can be found in the literature: Eulerian and Lagrangian. The former averages particle and fluid velocities * To whom correspondence should be addressed. E-mail: jose. [email protected]. † Castilla La Mancha University and Renewable Energy Research Institute. ‡ Chalmers University of Technology. § Carlos III University of Madrid.

over a volume larger than the particle size, and the resulting equations are similar to the equations of two interpenetrating fluids. In the latter, Newton’s equations of motion are applied for each individual particle, together with a model of energy interchange for the collisions between particles. The Eulerian approach needs some empirical correlations in order to close the problem,11 whereas in the Lagrangian approach, only systems with a limited amount of particles can be computed because of the exponential increase in the computational cost with the number of particles. Thus, further development is necessary to solve this problems, although CFD seems to be the most promising fluidized bed design tool for the future. This work tries to extend the works recently presented by Mu¨ller et al.12 and Laverman et al.,13 where PIV and digital image analysis were used to analyze the hydrodynamics of fluidized beds. In the present work, the particle velocity field is obtained using a PIV technique similar to the one used by Mu¨ller et al. and Laverman et al., but, in addition, the actual bubble geometry captured with a camera is exported to a finite element software, where the continuity and momentum equations for the gas flow can be solved numerically. In this way, the particle velocity obtained experimentally using PIV is used as input data in the numerical calculations. Figure 1 shows a scheme of the experimental-numerical approach of the paper. The gas flow equations proposed by Davidson’s model5,14 (gas continuity and Darcy law) are chosen to illustrate the application of the methodology proposed. The results obtained in a bubbling bed, with real bubble geometries and real particle velocity fields, are compared with the theoretical results obtained by the classical Davidson model (isolated circular bubbles and irrotational particle flow), showing qualitative agreement. In

10.1021/ie801720g CCC: $40.75  2009 American Chemical Society Published on Web 04/09/2009

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

5011

Figure 2. Characteristic curve of the distributor.

Figure 1. Scheme of the numerical-experimental approach followed in the paper.

addition, the equations solved numerically are modified to introduce non-Darcy effects in the calculations. Finally, the particle image velocimetry (PIV) results are also used to discuss the error committed in the numerical analysis when the voidage variation near bubbles is neglected, as it is the case of Davidson’s model. In the following, first the experimental setup and the PIV technique used to measure the particle velocity are briefly described. Thereafter, methodology and expressions used in the finite element computations are analyzed together with the results obtained. Then, the results obtained including non-Darcy effects in the computations are analyzed and compared with the previous results neglecting them. Finally, a study of the error committed when voidage variation around bubbles is neglected is presented and discussed together with the other results of the paper. Experimental Setup and PIV Measurements Experimental Setup. The experimental facility employed during the experiments was similar to the one described by Almendros-Iba´n˜ez et al.,15 that is, a 2D (110 × 60 × 0.5 cm3) fluidized bed. The bed was illuminated from the front while having a black background in order to get high contrast between particles and bubbles. The size of each frame captured with the high speed video camera was approximately 25 × 20 cm2, which allows a proper application of PIV technique to obtain the particle velocity field around a bubble with the given camera resolution of 1 megapixel. In all the pictures showed in this work, the level y ) 0 corresponds to the air distributor, and the height at which the bubble was captured is indicated in the vertical axis. The horizontal axis only indicates the scale. The pictures were captured from the central region of the bed with a rate of 250 fps and an exposure time of 1/5000 s. The particles were glass spheres with a mean diameter of dp ) 350 µm and a density of Fp ) 2500 kg/m3 (type B according to Geldart’s classification1). The height of the fixed bed was around 30-35 cm, and the excess gas velocity in all experiments was in the range U/Umf ≈ 1.1-1.3, thereby avoiding particles to be entrained out from the bed.

The distributor consisted of a perforated plate with 1 mm i.d. holes with a spacing of 1 cm, resulting in 1.57% open area. Figure 2 shows the ∆Pdist - U curve of the distributor. The superficial gas velocity was around 1 m/s, thus the distributor pressure drop was ∆Pdist ∼ 5000 Pa. During the experiments, the pressure drop of the gas crossing the bed was approximately ∆Pbed ) 4000 Pa; therefore, the bed and the distributor pressure drops were of the same order. Under these conditions, the plenum and the air supply system are uncoupled with the rest of the bed and the bubbles are small and uniformly distributed along the bed.3,16 PIV Measurements. The motion of bubbles and particles in 2D beds has been widely studied and discussed elsewhere.12,17-19 Raso et al.17 observed experimentally in a 2D bed that the movement of tracer particles were identical inside the bed and next to the wall. Link et al.19 studied numerically the motion of particles in a pseudo-2D geometry and only observed small differences between the first layer of solids and the whole bed movement. Therefore, the movement of the first layer of particles captured with the camera should be representative of the movement of all particles within the bed. In addition, the bubble diameter is much larger than the bed thickness; thus, any effects of the curvature in bubble thickness can be considered negligible.12 The particle velocity field around bubbles was measured applying the PIV technique on the emulsion phase. Note that cross-correlation is applied to the whole emulsion inside of the PIV window rather than to a group of dispersed particles, as is the case in traditional PIV applications. Thus, the velocity vector obtained for each window represents the mean velocity of the approximately 100 particles contained in the 16 × 16 pixels window. Three previous satisfactory studies using this technique in a bubbling bed can be found in the literature. Santana et al.20 measured the distribution of particle ejection velocity from erupting bubbles in a freely bubbling bed. Later on, Mu¨ller et al.12 extended the work of Santana et al.20 measuring the particle velocity distribution not only in the dome contour of the erupting bubble but also in the emulsion around the bubbles breaking at the bed surface. They concluded that the flow of particles is irrotational except in a small region close to the wake of the bubble. Both of them used a camera with a resolution of 512 × 480 ) 0.25 megapixel while a camera with a higher resolution of 1280 × 1024 ) 1.3 megapixel is used in the present study. With this resolution, a particle diameter to pixel size ratio of ∼2 is obtained in an area of ∼25 × 20 cm2. Mu¨ller et al. obtained a similar ratio but with a lower spatial resolution (∼15 × 15 cm2) and higher particle diameters. The pictures were

5012

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

Figure 3. (a) Original picture and (b) black and white picture obtained from a threshold value. The particle velocity is not computed in the white region of picture b as it is assumed to be free of particles. The scale is in centimeters.

obtained at a frequency of 250 fps and a shutter opening time of 200 µs in order to properly capture the bubble geometry and to avoid blurring by the motion of the bubbles. More recently, Laverman et al.13 combined PIV with digital image analysis to obtain time-averaged particle velocity profiles in 2D beds also using a 1.3 megapixel camera. The PIV software employed in this study is the same used by Mu¨ller et al.12 in their work, MATPIV 1.6.1.21 Also, the same iterative method reducing progressively the size window and the same filters (signal-to-noise, global, and local filters) used by them were applied in this work. The few dispersed particles inside of the bubble which rain down from the bubble roof have been neglected. Thus, the software has been modified in order to detect automatically the bubbles geometry and the freeboard using a threshold algorithm,22 and thereby, these regions were masked.13 Figure 3a shows the picture of one noninteracting slow bubble approximating the bed surface. By “non-interacting bubble” is meant a bubble in a freely bubbling bed located at least one bubble diameter away and with a gas recirculation which is not interacting with other bubbles. Thus, the bubble shown in Figure 3 is not an isolated bubble, but more bubbles are present in the bed at the same time although not captured within the 20 × 25 cm2 region framed. Figure 3b shows the result of applying a threshold value to the picture showed in Figure 3a. The bubble

Figure 4. (a) PIV results (blue vectors) and bubble velocity Ub ) 0.31 m/s (white vector) for the bubble shown in Figure 3a and b particle velocity magnitude contours. The scale is in centimeters, and velocities are in centimeters per second.

interior and the freeboard are detected and excluded from PIV analysis since these regions are considered free of particles in the computations. Mu¨ller et al.12 mixed white and black particles in order to obtain better contrast results. Nevertheless, the present study provides satisfactory results using only one kind of particles (around 90% of the vectors presented in each picture obtained a signal-to-noise ratio higher than 1.3) and no extra particles were needed. Laverman et al.13 also captured the particle velocity around bubbles with only one kind of particles. Figure 4a shows the PIV result obtained in the vicinity of the bubble shown in Figure 3a. The nonsymmetry of the result is due to interactions with other bubbles which do not appear in Figure 3a. These bubbles approximate the leading bubble from the bottom-right zone of the picture, and the gas flow from one bubble to the other drags the particles located in that zone, leading to higher particle velocity in this region. Also, the particles at the top of the bubble tend to move to right, as the bubble does, since the pressure gradient is higher in that region due to a lower level of the free surface of the bed. Figure 4b shows the magnitude of the particles velocity, where the higher velocity in the region close to the top of the bubble is distinguished more clearly. Thus, PIV measurements capture the real particle velocity field, taking into account the influence of other bubbles and the

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

5013

Table 1. Values and Uncertainties of Different Variables Computed Numerically 0.45 ( 2 × 10-3 5.5 × 10-10 ( 10-11 m2 0.45 ( 0.02 m/s 4000 ( 25 Pa

εmf ke Umf ∆Pbed

of the photographed region. In order to ensure this, some frames after the selected bubbles were observed in order to corroborate that there is enough space between bubbles and their influence can be neglected. A higher bed area can not be captured because of the limitation in spatial resolution for PIV measurements. In this paper, the equations proposed by the simplest Eulerian model (Davidson’s model5,14) for the gas flow were chosen to illustrate the application of the method presented. Davidson’s model has analytic solution using the equations of potential flow for circular, elliptical,9 or kidney shaped bubbles.8 In this work the bubble geometry is not artificial, as it is captured with a camera from a freely bubbling bed and the particle flow is not assumed irrotational but obtained experimentally from the PIV measurements. Therefore, the gas equations to solve are the following:

Figure 5. (a) Mesh created for the geometry shown in Figure 3a and (b) detail of the mesh around the bubble. The scale is in centimeters.

bed surface. This real particle velocity field is introduced in the finite element program, as it is explained in the following section. Numerical Analysis. With the experimental data obtained by means of videorecording and PIV technique, the gas velocity u and the pressure field pf were computed using a finite element software (Comsol Multiphysics23). First, the geometry of the bed is exported to this software, where an unstructured mesh is created automatically in the domain. This mesh is finer close to the boundaries (where a change in permeability k occurs) as Figure 5a and b shows. The mesh plotted in Figure 5a has ∼5 × 104 triangular elements. Figure 5a shows the addition of a rectangle, 15 cm in height, under the geometry captured by the camera, resulting in a computational bed geometry of the same height as the experimental bed geometry. This is done in order to enable proper implementation of the boundary condition at the bottom of the bed, as explained below. With this rectangle, the effect of other bubbles, which could be present in the bed under the region photographed, is not taken into account. Anyway, this effect is negligible if the other bubbles are not very close to the bottom

∇·u)0

(1)

k u˜ ) u - V ) - ∇ pf µ

(2)

where u˜ is the relative velocity between gas and particles, k is the permeability of the porous media, and pf is the fluid pressure. The voidage ε was assumed constant in the dense phase and equal to the one at minimum fluidization conditions εmf and bubbles were assumed free of particles (ε ) 1).13 The mass of particles was measured before they were introduced to the bed and a value of εmf ) 0.45 ( 0.02 was obtained for our experimental conditions. Note that eq 3 and the Ergun equation (eq 11, used later in the non-Darcy analysis) are valid when the voidage is close or equal to εmf. Close to the boundary of the bubble higher values of ε could be expected. A detailed study of the error resulting from this voidage variation is carried out in a following section. With respect to the pressure field, different expressions for the permeability through a porous media, k, can be found in literature. In this work the linear term of Ergun’s equation, also known as the Carman-Kozeny equation: ke )

εmf2(φdp)2 150(1 - εmf)2

(3)

has been used. For the experimental conditions, a value of ke ∼ 10-10 m2 was obtained. The regions of the domain free of particles (bubbles and freeboard) were modeled as regions with a permeability kb . ke. In this way, the gas pressure drop across such regions can be neglected and the pressure inside the bubbles was approximately constant. The measured voidage at minimum fluidization conditions has an uncertainty of (0.02 m/s. Using the law of propagation of uncertainty, the values and the uncertainties of the different variables are obtained and shown in Table 1. In order to study the influence of the permeability, kb, in the model, the simple case of an isolated circular bubble with varying bubble permeability, kb, and constant emulsion permeability, ke ) 10-10 m2, was analyzed. Figure 6a shows the nondimensional superficial gas velocity in an equatorial plane parallel to the distributor that divides the bubble into two equal parts, for different values of kb. As expected, the gas flow crossing the bubble (its center is situated at x ) 0) is higher

5014

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

Figure 7. Pressure distribution around one bubble according with Darcy’s law. The pressure is in Pascals, and the scale is in centimeters. The atmospheric pressure is assumed as reference and equal to 0.

Figure 6. (a) Nondimensional superficial gas velocity for several kb crossing one isolated circular bubble in its equatorial plane. (b) Detail of the flow in the interior of the bubble. In all cases, ke ) 10-10 m2, and the units of kb are squared meters.

than in the emulsion phase due to its higher permeability and increases asymptotically to a maximum value of U ) 2Umf, where Umf is the value far away from the bubble. Thus, the gas flow crossing the bubble is q ) 4UmfRb. This result is equal to the one obtained by Davidson and Harrison14 using the potential flow theory for the exchange between the bubble and the emulsion phase. As a consequence, the present numerical analysis agrees, in the simplest case of an isolated circular bubble, with the original results of Davidson’s model. Figure 6b shows a zoom of the flow in the region -Rb e x e Rb. As can be seen, no appreciable differences were observed for values of kb > 10-5 m2. Although it is not plotted, for high kb values the pressure inside of the bubble was constant and equal to the pressure in a point situated at the same height as the center of the bubble and far away from such a bubble. The boundary conditions chosen were atmospheric pressure in the freeboard and a pressure equal to the weight of particles at the bottom of the bed. The same boundary conditions were assumed by Croxford2 in order to obtain the pressure field in a fluidized bed. For the lateral limits, a boundary condition of no-penetration of gas, i.e. the pressure gradient perpendicular to the boundary is equal to zero, was imposed. Nevertheless, Croxford showed that this lateral boundary condition has little effect on the rest of the flow. Note that by assuming the pressure at the bottom of the bed equal to the hydrostatic pressure, the effect of the gas discharge in the distributor is neglected, i.e. a perfect gas flow distribution is

assumed. This boundary condition is valid when the distributor pressure drop is high enough and the bed and the air-supply system are uncoupled, as was explained in the Experimental Setup section. In addition, the flow at the region close to the bottom of the bed was visually inspected, observing an apparently homogeneous distribution of the flow. If the pressure drop across the distributor was not high enough, a different behavior could be observed in the bed, as the single or exploding bubble regimes observed by Johnsson et al.3 In this situation, the assumption of constant pressure at the bottom of the bed could not be correct (i.e., the static pressure as well as the gas flow will fluctuate). From the numerical calculations (using the PIV measurements as input data in the numerical analysis), the gas velocity u and the pressure field pf were obtained solving eqs 1 and 2. The following show the numerical results obtained for slow and fast bubbles, as well as erupting and interacting ones. Slow and Fast Bubbles. The pressure distribution around a slow bubble obtained solving eqs 1 and 2 is shown in Figure 7. Distortion of the isobars by the presence of the bubble can be seen, whereas the pressure in its interior is constant. Also the fact that the bed surface is not horizontal forces the isobars to adapt to it. The pressure drop in the freeboard is negligible respect to the pressure drop across the bed due to the differences in the permeabilities of both media. The interstitial gas velocity vectors and their streamlines are plotted in Figure 8a and b, respectively, as they were observed by a noninertial viewer. Deviation of the gas streamlines trajectories in the vicinity of the bubble due to the more favorable path across the region of higher permeability is observed. Nevertheless, more interesting is to plot the velocity field, which is viewed by an observer moving with the bubble, as shown in Figure 8c and d. These figures show, for this slow bubble (Ub < umf), two small regions by both sides of the bubble where the gas recirculates and how the rest of the gas traverses the bubble from bottom to top, in agrement with the flow predicted by Davidson’s model in the slow bubble case. The radius of the percolation circle that defines the limits of the lateral recirculation zones can be obtained using the 2D Davidson’s model as the following: Rp ) Rb

U +u | U  - u | ) 3.12 cm b

mf

b

mf

(4)

where Rb has been calculated as the radius of the equivalent circular bubble with the same area. The circular bubble model seems to overestimate Rp which has a value of approximately ∼2.5 cm. This difference is owing to the noncircular geometry of the bubble and the nonaxissymmetric particle velocity field.

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

A different behavior is predicted by Davidson’s model when the bubble velocity is higher than the interstitial gas velocity. These fast bubbles are characterized by a cloud of gas recirculating in the bubble without penetration of external gas inside of this recirculation region. Figure 9a and b shows the velocity vectors and the streamlines viewed by a stationary observer, and Figures 9c and d shows that viewed by an observer moving with the bubble for a noninteracting fast bubble with a ratio between the bubble velocity and the interstitial velocity of Ub/umf ) 1.3. Typical recirculation regions can be observed in both sides of the bubble, although an important part of the gas crosses the bubble from bottom to top and is not recirculated. This fact is the result of the nonvertical bubble ascent (note the appreciable horizontal component of the bubble velocity represented with a white colored vector in Figure 9a and c). In fact, Davidson’s model shows a similar behavior for a circular bubble if an horizontal component is added to the bubble velocity. It is well-known from Davidson’s theory that the gas stream function for a circular bubble can be obtained as follows:8

5015

ψg ) ψg0 + ψp where

(5)

( )

ψg0 ) -umf 1 +

Rb2

r sin(θ) (6) r2 is the gas stream function with Ub ) 0, i.e., ψg0 represents the percolation of fluid through a fixed bed into a circle of radius Rb, and it is not affected by the bubble velocity. In contrast, the particle stream function ψp must be modified in order to take into account that the bubble is moving with a velocity forming and angle R with the vertical, that is

( )

ψp ) Ub 1 -

Rb2

r sin(θ + R) r2 Combining eqs 5-7 leads to following expression

[ ]

ψg ) (Ub cos(R) - umf)r sin(θ) 1 -

j2 A + r2

(7)

[ ]

Ub sin(R)r cos(θ) 1 -

Rb2 r2

(8)

Figure 8. (a) Gas velocity vectors u and (b) streamlines as viewed by a stationary observer and (c and d) as viewed by an observer moving with the bubble. Slow bubble with Ub/umf ) 0.69 and Ub ) 0.31 m/s. The scale is in centimeters.

5016

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

Figure 9. (a) Gas velocity vectors u and (b) streamlines as viewed by a stationary observer and (c and d) as viewed by an observer moving with the bubble. Fast bubble with Ub/umf ) 1.3. White vector indicates the bubble velocity Ub ) 0.58 m/s. The scale is in centimeters.

j 2 ) [(Ub cos(R) + umf)/(Ub cos(R) - umf)]Rb2 and the where A radial and tangential gas velocities can be obtained as follows: ur )

uθ )

cos(θ + R) -1 ∂ψg ) cos(θ) Ub + umf × r ∂θ cos(θ) Rb2 Ub cos(θ + R) - umf cos(θ) (9) Ub cos(θ + R) + umf cos(θ) r2

(

)

[

]

∂ψg sin(θ + R) ) sin(θ) Ub + umf × ∂r sin(θ) 2 Ub sin(θ + R) - umf sin(θ) Rb + (10) 2 U r b sin(θ + R) + umf sin(θ)

(

[

)

]

The first term of the right side of eq 8 represents the streamlines due to the vertical component of the bubble velocity and has the same form as the original Davidson’s equation changing Ub by Ub cos(R). The second term represents the irrotational flow of particles around a circular bubble moving with an horizontal velocity Ub sin(R). Note that eqs 8-10 recover the original Davidson’s expressions for R ) 0. Figure 10 shows the streamlines obtained from the stream function 8 for a circular bubble with a velocity equal to the velocity of the bubble shown in Figure 9. A qualitative agreement between the simple case of a circular bubble with the one shown in Figure

Figure 10. Gas streamlines according with eq 8 for a fast bubble with Ub ) 0.58 m/s, umf ) 0.446, and R ) 21.8 °.

9d is observed: two recirculation regions at both sides of the bubble but also part of the gas that traverses and leaves the fast bubble. However, not only does the direction of the bubble velocity affect the gas velocity field, but also the bubble geometry significantly influences the gas flow. This can be seen in Figure 9 where higher velocities are obtained in the left-bottom zone of the bubble because of the favorable path created by the bubble

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

5017

Figure 11. Gas velocity vectors u and streamlines for an erupting bubble viewed by a stationary observer. The scale is in centimeters.

geometry in that region. This fact can not be observed in a circular bubble. Thus, the throughflow crossing the bubble is also influenced by the bubble geometry. Erupting and Interacting Bubbles. Previous studies about erupting bubbles can be found in the literature focused either on the velocity of ejected particles12,15,20 or on the dome evolution and the gas flow through the bubble.24-28 Although the particle velocity field was measured experimentally in bubbling beds using PIV, the few studies available about the throughflow in erupting bubbles are limited to isolated circular

bubbles or artificial cavities at the bed surface rather than in real geometries. Figure 11 shows the evolution of the velocity field and the streamlines viewed by an stationary observer for one bubble eruption at the bed surface. At t ) 0 the bubble is approximating to the bed surface forming the typical dome. At t ) 40 ms the dome breaks at a certain point due to the irregularities which appear in the dome when its thickness decreases and the stalactites funnel the particles situated in the external surface of the dome.12,28 In this instant, most of the throughflow leaves

5018

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

Figure 12. (a-c) Streamlines viewed by a stationary observer. (d-f) Streamlines viewed by an observer moving with the leading bubble. The scale is in centimeters.

the bubble through the small aperture opened at the breaking dome. The magnitude of the gas velocity is very high in this region, and the rest of the dome decelerates because the gas that was dragging it is drastically diminished and dome particles are moving only because of their inertia. Finally at t ) 68 ms the dome collapses and the gas velocity decreases due to the higher section opened for the gas flow at the dome aperture. In addition, during the dome evolution shown in Figure 11, a small bubble is observed at the bottom-right region of the pictures. This bubble is seen to follow the path opened by the leading one due to the more favorable pressure gradient, although they do not coalesce before the leading bubble breaks the bed surface. Figure 12 shows the evolution of one erupting bubble while another bubble is coalescing below the leading one. At t ) 0 the bubble is isolated and it approximates, ascending vertically, to the bed surface with a velocity Ub > umf, as it is deduced from the gas recirculation vortices at both sides of the bubble observed in Figure 12d, which are typical of fast bubbles. At t ) 92 ms the leading bubble continues ascending through the bed with a velocity Ub > umf, but now another elongated bubble is approximating to the leading bubble. The more favorable pressure gradient through the path opened by the first bubble provokes an acceleration in the trailing bubble, which reaches a higher velocity. In this situation, a single gas cloud encompassing both bubbles can be seen, in agreement with flow patterns of gas around a pair of bubbles aligned vertically observed by Shichi et al.29 (also reproduced in the first edition

of Kunii and Levenspiel30). In addition, the peculiar geometry of both bubbles provoques additional smaller gas recirculations in the right side of both bubbles. Between Figure 12a and b, the throughflow crossing the leader bubble ranges between 2 and 2.6 times Umf. Finally, at t ) 124 ms, the leading bubble decelerates and becomes a slow bubble (Ub < umf). The gas recirculations disappears, and the roof of the bubble crumbles because of the rain of particles from the bubble roof in form of stalactites. The sequence of pictures of Figure 12 also shows the change in the geometry of both bubbles. The leading bubble, initially slightly elongated, begins to grow and to expand in the horizontal direction as it breaks the bed surface. This bubble corresponds to the collapsed dome bubbles observed by Almendros-Iba´n˜ez et al.15 This kind of bubble is characterized by a very low bubble velocity and nonprojected dome particles at eruption; rather, dome particles rain off and return to the dense phase. They usually appear when almost the entire bubble is over the bed surface at the eruption instant. In this situation, throughflow crossing the bubble and, consequently, the particle ejection velocity are low.24 During the last steps of the erupting process, the gas flow crossing the leader bubble is reduced around 1.5Umf. The trailing bubble, which coalesces below the leading one, also changes its geometry: it initially elongates due to the acceleration during the approximation to the wake of the leading bubble. This picture is typical of coalescing bubbles. Later on, as the bubble decelerates, its geometry changes to an ap-

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

5019

31-33

beds have demonstrated that Darcy’s law is valid only in the creeping flow regime, i.e. at low Reynolds numbers. Ergun’s equation34 (1 - εmf)2 µU 1 - εmf FfU2 ∆pf ) 150 + 1.75 L ε 3 (φd )2 ε 3 φdp mf

Figure 13. Graphical representation of the nondimensionalized Ergun’s equation. The dashed line represents eq 12, and the solid line represents only the lineal term of eq 12. The cross represents the operating conditions at minimum fluidization conditions, and the vertical lines indicate the range of variation of the local Reynolds number for the experimental conditions of this work.

proximately round bubble, due to the increase of the wake section of the leading bubble. Numerical Analysis Including Non-Darcy Effects. Inertial (or non-Darcy) effects are not taken into account when Darcy’s law (eq 2) is assumed as the gas momentum equation. Numerous studies about the flow through porous media and packed

p

(11)

mf

combines a linear term (Carman-Kozeny term), which predicts properly the pressure drop when the viscous forces dominate, and a quadratic term (Forchheimer term) which is valid when the inertial forces are more important than the viscous ones. It is generally accepted, and has also been experimentally demonstrated, that the sum of both terms predicts properly the pressure drop when both forces are of the same order of magnitude, i.e., at intermediate Reynolds numbers.35 Equation 11 can be rearranged and nondimensionalized according to

( ) ∆pf L

*

)

85.714 +1 Re

(12)

where (∆pf/L)* ) (∆P/L)[(εmf3dpφ)/(1.75(1 - εmf)FfU2)] is a nondimensional pressure gradient and Re is the Reynolds number defined by FfU Re )

φdp (1 - εmf) µ

(13)

Figure 14. (a) Gas streamlines of one slow bubble (Ub/umf ) 0.69) viewed by an observer moving with the bubble including non-Darcy effects and (b) relative difference in the magnitude of u. The scale is in centimeters.

Figure 15. (a) Gas streamlines of one fast bubble (Ub/umf ) 1.3) viewed by an observer moving with the bubble including non-Darcy effects and (b) relative difference in the magnitude of u. The scale is in centimeters.

5020

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

Figure 16. (a) Gas streamlines of one erupting bubble, in the instant when the bubble breaks the bed surface, viewed by a stationary observer including non-Darcy effects and (b) relative difference in the magnitude of u. The scale is in centimeters.

Equation 12 is plotted in Figure 13, which shows that approximately around Re ∼ 1-10 Ergun’s equation departs from the lineal relation proposed by Darcy’s law. For Re g 104, the value of (∆pf/L)* is constant and independent of the Reynolds number. In this figure, the cross indicates the Reynolds number obtained at minimum fluidization conditions: Remf ) 9. Nevertheless, this is a mean value and the gas velocity can vary notably in the vicinity of the bubbles. For the experimental conditions of this work, the local Reynolds number ranges between 4 and 40. In regions of higher Reynolds number, the influence of non-Darcy effects will be significant. In order to take into account non-Darcy deviations, Ergun’s equation can be used instead of Darcy’s law. In this way, the magnitude of the pressure gradient can be obtained from the following equation | ∇ pf| ) Bu˜ + Au˜2

(14)

where A ) 1.75{(Ff/dpφ)[(1 - εmf)/εmf]} and B ) 150 {(µ/dp2φ2)[(1 - εmf)2/εmf2]} are constants obtained from eq 11. The magnitude of the relative gas velocity can be obtained from (14) as follows -B + √B2 + 4A| ∇ pf| (15) 2A and its components are proportional to the pressure gradient in each direction according to u˜ )

∂pf ∂x u˜x ) u˜ ; | ∇ pf| -

∂pf ∂y u˜y ) u˜ | ∇ pf| -

(16)

Figure 14a shows the streamlines for the slow bubble case shown in Figure 8 obtained taking into account non-Darcy effects. Both streamlines plots look quite similar, and the nonDarcy effects do not seem to affect appreciably the gas path, although some differences are observed in the magnitude of the gas velocity. The influence of the non-Darcy effects on the magnitude of the gas velocity is analyzed through a plot of the parameter ∆u, which is defined by eq 17 and represents the relative difference between the magnitude of the gas velocity including or not non-Darcy effects. ∆u )

uD - uD-F uD

(17)

Figure 14b shows the results obtained for the slow bubble case. The maximum differences (∆u ) 0.12) are observed in the nose and in the wake of the bubble, where the magnitude of u is higher (and thereby the local Reynolds number) due to the preferential path for the gas flow created by the bubble. In contrast, differences are small at both sides of the bubble. Figure 15 shows the fast bubble case. Similar results to the ones obtained in the slow bubble are observed when non-Darcy effects are included: higher differences ∆u in the region close to the top of the bubble and lower in the recirculation regions. Although in this case, the maximum difference (∆u ∼ 0.3) is higher, due to the higher bubble velocity. In both slow and fast bubbles, as expected, the values of uD are higher than those of uD-F. For the same pressure gradient, the magnitude of the gas velocity is higher neglecting Forchheimer effects. The difference increases with the Reynolds number, although the streamlines in the vicinity of the bubble do not seem to be affected appreciably. Figure 16 shows the erupting bubble case at the instant when the bubble breaks the bed surface. Most of the gas crossing the bubble reaches the freeboard through the small aperture opened on the top of the dome contour. Although the velocity in this aperture is very high, the value of ∆u in this region is not high because the interior of the bubble is connected to the freeboard and both are modeled as regions of very high permeability and the nonlinear effects are not important in these regions. The streamlines and gas velocity differences ∆u taking into account non-Darcy effects for two interacting bubbles at the instant when the streamlines of the leading bubble interact with the second bubble can be seen in Figure 17. As in the previous cases, highest differences in the gas velocity are observed in the top of the leading bubble but also in the region situated between the bubbles. In summary, under the experimental conditions of our work (Reynolds number ranging between 4 and 40), inertial effects have only influence in the magnitude of the gas velocity in regions where the gas velocity, and consequently the local Reynolds number, are higher. These regions are typically the wake and the nose of the bubble. Nevertheless, the gas streamlines are not appreciably affected by non-Darcy effects, and the general gas path does not change.

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

5021

Figure 17. (a) Gas streamlines of two interacting bubbles viewed by an observer moving with the leading bubble including non-Darcy effects and (b) relative difference in the magnitude of u. The scale is in centimeters.

Influence of Voidage Variation around Bubbles and Discussion. Influence of Voidage Variation around Bubbles. The numerical results given in the previous sections were obtained assuming a constant value for the voidage εmf in the dense phase. Nevertheless, several studies in the literature36-41 showed that there is small region around the bubbles where the voidage decreases from the value at the bubble contour to a value close or equal to εmf. The particle velocity field measured using PIV can be used to analyze the influence of the voidage distribution in the final results. In this way, the general continuity equation for the solid phase is6,7 ∂ ((1 - ε)) + ∇ · ((1 - ε)V) ) 0 ∂t

(18)

which can be rearranged as follows -1 Dε + ∇ ·V)0 (1 - ε) Dt

(19)

where Dε/Dt ) ∂ε/∂t + V · ∇ε represents the substantive derivative of ε. The substantial derivative of the voidage should only be significant in the region adjacent to the bubbles, whereas it should be negligible in the rest of the bed, where ε is approximately constant. Therefore, far away from the bubble, the continuity equation reduces to its classical form for incompressible flow: ∇·V)0 (20) The difference between eq 20 and the actual particle velocity field obtained using PIV shows the regions of the bed where significant voidage gradients exist. Figure 18 shows the result obtained calculating the divergence of the particle velocity ∇ · V from the PIV measurements for the slow bubble case shown in Figures 4 and 8. As expected, ∇ · V = 0 far from the bubble, but in the region adjacent to the bubble, ∇ · V * 0 due to the importance of the porosity gradients (first term of eq 19). Therefore, the voidage influence is important only in a region of a few millimeters around the bubble. Similar results were obtained for the rest of the bubbles shown in this work.

Figure 18. Divergence of the particle velocity field ∇ · V obtained from PIV measurements. The scale is in Hertz.

Figure 19 shows the divergence of the particle velocity in a small region at the nose of the erupting bubble of Figure 11. Higher voidage values are expected in this region, and consequently, higher differences in the divergence than in the bubble showed in Figure 18 are observed. Note the difference in the scales between Figures 18 and 19, which expresses the higher voidage gradients at the bubble nose of the erupting bubble. In addition, the highest values of the divergence of the particle velocity are at the top of the bubble nose (the region marked with a square in Figure 19), which is the point where the throughflow crossing the bubble breaks the bed surface 40 ms later (see Figure 11c) due to the lower resistance to the gas flow. In another work Almendros-Iba´n˜ez et al.41 studied in more detail the influence of the voidage variation around bubbles on the flow patterns around the bubble under similar experimental conditions, although using a different experimental technique. This consists of illuminating with diffuse light from the rear of the bed and transforming the gray levels around the bubble into voidage values. They observed that, although the throughflow crossing the bubble increases due to the higher permeability

5022

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

of bubbles traverses the bed as throughflow crossing the bubbles. They also concluded that the throughflow is not dependent on the bubble rise velocity or particle velocities but on the bubble shape and bubble distribution. Here, the instantaneous throughflow for each frame (rather than the time-averaged throughflow, which is the data usually presented in literature43,44) is computed for real bubble geometries. In agreement with Valenzuela and Glicksman,42 bubble geometry is shown to play an important role in the throughflow, especially in the corners of the bubble contour (see Figure 9c) where high gas velocities are observed. The numerical results of Almendros-Iba´n˜ez et al.41 corroborate the importance of the bubble geometry. Conclusions

Figure 19. Detail of the divergence of the particle velocity field ∇ · V obtained from PIV measurements at the bubble nose of the erupting bubble shown in Figure 11a. The scale is in Hertz.

region around the bubble, the general gas flow pattern does not change appreciably. Discussion The results shown in this work were obtained assuming a constant voidage in the dense phase ε ) εmf and ε ) 1 within the bubbles. Different techniques can be found in the literature to measure voidage around bubbles. For example, Lockett and Harrison36 and Nguyen et al.37 used capacitance probes in 2D beds and Yates et al.39 used X-rays to obtain voidage values in a 3D bed with injected bubbles. None of these techniques has enough spatial resolution to accurately measure voidage in the region adjacent to the bubble where the voidage gradients are important (it follows from our results that a spatial resolution of a few millimeters is required). In a previous work, Almendros-Iba´n˜ez et al.41 obtained voidage distribution by illuminating a 2D bed from the rear, which is incompatible with simultaneous PIV measurements. Nevertheless, although there is no technique at present that permits us to obtain particle velocities and voidage values simultaneously with enough spatial resolution, AlmendrosIba´n˜ez et al.41 showed that the influence of the voidage field on the gas flow pattern is small. In the same way, Laverman13 also assumed ε ) 1 within the bubble and ε ) εmf in the dense phase in their hydrodynamic study using PIV to obtain timeaveraged solid mass flows. The methodology shown in this work can also be used to validate the results obtained by CFD simulations. In this way, the PIV results can be used not only to compare the particle velocity field obtained numerically, but also as input to obtain a gas velocity field. This gas velocity field can be used to obtain numerically a particle velocity field, which should be compared with the PIV results. Thus, the PIV results can be used to check the numerical results of both particle and gas phases. In the present case, this is not posible due to the fact that the particle velocity field is independent of the gas velocity in the Davidson model, i.e. only depends on the geometry and bubble velocity. Valenzuela and Glicksman42 computed the throughflow crossing bubbles for an array of spherical bubbles distributed in a bed. They showed that, if the porosity in the emulsion phase remains at the value corresponding to minimum fluidization conditions, all excess gas over Umf which is not visible in form

The combination of digital image analysis, PIV, and FEM makes it possible to obtain the velocity fields of the particle phase (experimentally) as well as the gas phase (through numerical simulations) around bubbles in fluidized beds. This method allows evaluation of different models for the gas phase flow and opens for future experimental validation of simulated pressure fields. The solid-gas flow has been analyzed around different types of bubbles: slow, fast, erupting, and interacting. It is also shown that the equations of the classical Davidson model (developed from potential flow theory for circular bubble geometries) can be applied to real, noncircular bubbles with experimentally determined particle velocities. In addition, the effect of the nonvertical ascent of fast bubbles is analyzed with respect to how these bubbles interchange gas with the main stream, and the influence of the voidage around bubbles on the final results is discussed. Acknowledgment This work has been supported by the National Energy Programme of the Spanish Ministry of Education under the project number ENE2006-01401. Nomenclature A ) constant defined in eq 14 [kg/m4] j ) constant defined in eq 8 [m2] A B ) constant defined in eq 14 [kg/(s · m3)] dp ) particle diameter [m] k ) permeability of the medium [m2] kb ) permeability of the bubble phase [m2] ke ) permeability of the dense phase (emulsion) [m2] L ) length of the bed [m] pf ) fluid pressure [Pa] q ) gas flow through a bubble [m2/s] Rb ) bubble radius [m] Rp ) radius of the percolation circle around bubbles, defined by eq 4 [m] Re ) Reynolds number defined by eq 13 [-] r ) radial coordinate t ) time [s] U ) superficial gas velocity [m/s] Ub ) bubble velocity [m/s] Umf ) superficial gas velocity at minimum fluidization conditions [m/s] u ) interstitial gas velocity [m/s] umf ) interstitial gas velocity at minimum fluidization conditions [m/s] ur ) radial gas velocity [m/s] uθ ) tangential gas velocity [m/s] u˜ ) relative gas velocity defined in eq 2 [m/s]

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009 V ) particle velocity [m/s] x ) horizontal coordinate y ) vertical coordinate R ) angle formed by the bubble velocity and the vertical [rad] ∆Pdist ) pressure drop of the gas crossing the distributor [Pa] ∆pf ) fluid pressure drop [Pa] ∆u ) nondimensional gas velocity difference defined by eq 17 [-] ε ) void fraction [-] εmf ) void fraction at minimum fluidization conditions [-] θ ) tangential coordinate µ ) fluid dynamic viscosity [Pa · s] Ff ) fluid density [kg/m3] Fp ) particle density [kg/m3] φ ) particle sphericity [-] ψg ) gas stream function ψg0 ) gas stream function with Ub ) 0 ψp ) particle stream function ()D ) calculated using the Darcy equation ()D-F ) calculated using the Darcy-Forchheimer equation ()x ) x component ()y ) y component

Literature Cited (1) Geldart, D. Types of gas fluidization. Powder Technol. 1973, 7, 285. (2) Croxford, A. J. Investigation of dynamics and control of fluidised beds using pressure measurements. PhD thesis, University of Bristol, 2006. (3) Johnsson, F.; Zijerveld, R. C.; Schouten, J. C.; van den Bleek, C. M.; Leckner, B. Characterization of fuidization regimes by time-series analysis of pressure fuctuations. Int. J. Multiphase Flow 2000, 26, 663. (4) Pallare`s, D.; Johnsson, F. A novel technique for particle tracking in cold 2-dimensional fluidized beds-simulating fuel dispersion. Chem. Eng. Sci. 2006, 61, 2710. (5) Davidson, J. F. Symoposium on Fluidisation-Discussion. Trans. Inst. Chem. Eng. 1961, 39, 230. (6) Jackson, J. R. The mechanics of fluidised beds: Part II: The motion of fully developed bubbles. Trans. Inst. Chem. Eng. 1963, 41, 22. (7) Murray, J. D. On the mathematics of fluidization. Part 2. Steady motion of fully developed bubbles. J. Fluid Mech. 1965, 22, 57. (8) Collins, R. An extension of Davidson’s theory of bubbles in fluidized beds. Chem. Eng. Sci. 1965, 20, 747. (9) Clift, R.; Grace, J. R.; Cheung, L.; Do, T. H. Gas and solids motion around deformed and interacting bubbles in fluidized beds. J. Fluid Mech. 1972, 51, 187. (10) van Wachem, B. G. M.; Almstedt, A. E. Methods for multiphase computational fluid dynamics. Chem. Eng. J. 2003, 96, 81. (11) Jackson, J. R. The dynamic of fluidized particles; Cambridge University Press: Cambridge, 2000. (12) Mu¨ller, C. R.; Davidson, J. F.; Dennis, J. S.; Hayhurst, A. N. A study of the motion and eruption of a bubble at the surface of a two-dimensional fluidized bed using Particle Image Velocimetry (PIV). Ind. Eng. Chem. Res. 2007, 46, 1642. (13) Laverman, J. A.; Roghair, I.; van Sin Annaland, M.; Kuipers, H. Investigation into the hydrodynamics of gas-solid fluidized beds using particle image velocimetry coupled with digital image analysis. Can. J. Chem. Eng. 2008, 86, 523. (14) Davidson, J. F.; Harrison, D. Fluidised particles; Cambridge University Press: Cambridge, 1963. (15) Almendros-Iba´n˜ez, J. A.; Sobrino, C.; de Vega, M.; Santana, D. A new model for ejected particle velocity from erupting bubbles in 2D fluidized beds. Chem. Eng. Sci. 2006, 61, 5981. (16) Sasic, S.; Johnsson, F.; Leckner, B. Fluctuations and waves in fluidized bed systems: The influence of the air-supply system. Powder Technol. 2005, 153, 176. (17) Raso, G. G.; Volpicelli, G.; Maitz, C. Sul movimento do bolle di gas in letti di particelle fluidizzati in colonne rettangolari. Quaderni dell’ingegnere chimico italiano 1965, 1, 157. (18) Gabor, J. B. Wall effects on fluidized particle movement in a twodimensional column. In International Symposium on Fluidization; Drinkenburg, A. A. H., Ed.; Netherlands University Press: Amsterdam, 1967; Eidenhoven, June 6–7, 1967; p 230.

5023

(19) Link, J.; Zeilstra, C.; Niels, D.; Kuipers, H. Validation of a discrete particle model in a 2D spout-fluid bed using non-intrusive optical measuring techniques. Can. J. Chem. Eng. 2004, 82, 30. (20) Santana, D.; Nauri, S.; Acosta, A.; Garcı´a, N.; Macı´as-Machı´n, A. Initial particle velocity spatial distribution from 2D erupting bubbles in fluidized beds. Powder Technol. 2005, 150, 1. (21) Sveen, J.P. http://www.math.uio.no/∼jks/matpiv (last modified in August 2004; accessed in 2007). (22) Otsu, N. A threshold selection method from gray-level histograms. IEEE Trans. Syst. Man Cybern. 1979, 9, 62. (23) www.comsol.com (accessed Sep 2008). (24) Glicksman, L. R.; Yule, T. W. Gas throughflow in a bubbling fluidized bed. In Proceedings of the Engineering Foundation Conference, Fluidization V, Elsinore (Denmark), May 18–23, 1986; Bear, J., Corapcioglu, M. Y., Eds.; p 103. (25) Glicksman, L. R.; Yule, T. W. Prediction of the particle flow conditions in the freeboard of a freely bubbling fluidized bed. Chem. Eng. Sci. 1995, 50, 69. (26) Levy, E. K.; Chen, H. K.; Radcliff, R.; Caram, H. S. Analysis of gas flow through erupting bubbles in a gas-fluidized bed. Powder Technol. 1988, 54, 45. (27) Gera, D.; Gautam, M. Variation of throughflow velocity in a 2D rising bubble. Powder Technol. 1994, 79, 257. (28) Almendros-Iba´n˜ez, J. A.; Sa´nchez-Delgado, S.; Sobrino, C.; Santana, D. Experimental observations on the different mechanism for solid ejection in gas-solid fluidized beds. Chem. Eng. Process. 2009, 48, 734. (29) Shichi, R.; Mori, S.; Muchi, I. Interaction between two bubbles in gaseous fluidization. Kagaku Kogaku 1968, 32, 343. (30) Kunii, D.; Levenspiel, O. Fluidization Engineering; John Wiley & Sons: New York, 1969. (31) Dybbs, A.; Edwards, R. V. A new look at porous media fluid mechanics - Darcy to turbulence. In Fundamentals of Transport Phenomena in Porous Media; Bear, J., Corapcioglu, M. Y., Eds.; Martinus Nijhoff Publishers: Dordrecht, 1984; p 199. (32) Bear, J. Dynamics of fluids in porous media; American Elsevier: New York, 1972. (33) Zenz, Z.; Grigg, G. A Criterion for Non-Darcy Flow in Porous Media. Transp. Porous Media 2006, 63, 57. (34) Ergun, S. Fluid flow through packed columns. Chem. Eng. Prog. 1952, 48, 89. (35) Bird, R. B.; Stewart, W. E.; Lightfoot, E. N. Transport phenomena, 2nd ed.; John Wiley & Sons, Inc.: New York, 2002. (36) Lockett, M. J.; Harrison, D. The distribution of voidage fraction near bubbles rising in gas-fluidized beds. In International Symposium on Fluidization; Drinkenburg, A. A. H., Ed.; Netherlands University Press: Amsterdam, 1967; Eidenhoven, June 6–7, 1967; p 257. (37) Nguyen, X. T.; Leung, L. S.; Weiland, R. H. On void fractions around a bubble in a two dimensional fluidized bed. In Proceedings of the International Congress on Fluidization and its Applications, Toulouse, October 1–5, 1973; Angelino, H., Couderc, J. P., Gilbert, H., Lagarie, C., Eds.; p 230. (38) Collins, R. A model for the effects of the voidage distribution around a fluidization bubble. Chem. Eng. Sci. 1989, 44, 1481. (39) Yates, J. G.; Cheesman, D. J.; Sergeev, Y. A. Experimental observations of voidage distribution around bubbles in a fluidized bed. Chem. Eng. Sci. 1994, 49, 1885. (40) Buyevich, Y. A.; Yates, J. G.; Cheesman, D. J.; Wu, K. T. A model for the distribution of voidage around bubbles in a fluidized bed. Chem. Eng. Sci. 1995, 50, 3155. (41) Almendros-Iba´n˜ez, J. A.; Pallare`s, D.; Johnsson, F.; Santana, D. Voidage distribution around bubbles in a fluidized bed: influence on throughflow. Powder Technol., submitted for publication. (42) Valenzuela, J. A.; Glicksman, L. R. Gas flow distribution in a bubbling fluidized bed. Powder Technol. 1985, 44, 103. (43) Grace, J. R.; Clift, R. On the two-phase theory of fluidization. Chem. Eng. Sci. 1974, 29, 327. (44) Johnsson, F.; Andersson, S.; Leckner, B. Expansion of a freely bubbling fluidized bed. Powder Technol. 1991, 68, 117.

ReceiVed for reView November 11, 2008 ReVised manuscript receiVed March 1, 2009 Accepted March 17, 2009 IE801720G