Assessment of Model Formulations in the Discrete Particle Simulation

8378

Ind. Eng. Chem. Res. 2004, 43, 8378-8390

Assessment of Model Formulations in the Discrete Particle Simulation of Gas-Solid Flow Y. Q. Feng and A. B. Yu* Centre for Simulation and Modeling of Particulate Systems, School of Materials Science and Engineering, The University of New South Wales, Sydney NSW 2052, Australia

Discrete particle simulation has been recognized as a useful numerical technique for elucidating the fundamentals of granular matter. For gas-solid two-phase flow in fluidization, such simulations are achieved by combining the discrete flow of the particle phase with the continuum flow of the gas phase. However, differences exist in the actual implementation of this idea in the literature. This paper attempts to rationalize this matter by discussing important aspects including the governing equations in relation to the so-called models A and B, which use different treatments of pressure drop in the well-established two-fluid model, different coupling schemes between the gas and solid phases, and different equations for quantifying the particle-fluid interaction. For the purpose of quantitative analysis, gas fluidization of binary mixtures of particles is simulated with different model formulations, and a comparison of the results in terms of flow pattern and mixing/segregation kinetics shows a significant difference. Physical experiments are then conducted under similar conditions to assess the two model formulations, and the results suggest that the model B treatment is favored. The reason for the difference is also discussed in terms of the particle-fluid interaction force. On this basis, a method to reduce the difference between the two model formulations is proposed and tested. Introduction Gas-solid flow in fluidization can be found in many industries and has been the subject of intensive research in the past. Numerical simulation has played an important role in the development of a comprehensive understanding of the fundamentals of this phenomenon. Generally speaking, such simulations can be done by either continuum or discrete methods, depending on the time and length scales employed for solid phase. The former is often based on the so-called two-fluid model (TFM) (see ref 1, for example), whereas the latter is usually based on the discrete element method (DEM)2 or, more generally, the discrete particle method (DPM). The discrete approach, when coupled with the gas flow, which is often solved by computational fluid dynamics (CFD), gives the so-called DPM-CFD or combined continuum and discrete model (CCDM), representing the fact that the nature of such modeling is a combination of continuum gas and discrete solid flows.3-5 In the CCDM, the motion of discrete particles is obtained by solving Newton’s second law of motion as used in the DEM and the flow of continuum gas by solving the Navier-Stokes equations based on the concept of local average as used in the TFM, with their coupling through a particle-fluid interaction force. The CCDM can generate detailed particle-scale information, such as the trajectories of and forces acting on individual particles, which is key to elucidating the mechanisms governing the complicated flow behavior. Such information cannot be obtained in detail with the TFM or with the current experimental techniques. Therefore, it is particularly suitable for fundamental research. With the rapid development of computer technology, the CCDM has been increasingly used by various investigators to * To whom correspondence should be addressed. Tel.: 612-9385-4429. Fax: 61-2-9385-5956. E-mail: [email protected].

Table 1. Model Formulations, Schemes for Coupling between Particle and Fluid, and Equations for Calculating the Particle-Fluid Interaction Force Available in the Literature governing equations two-way coupling

scheme 1 scheme 2 scheme 3 particleErgun40 and fluid Wen and Yu41 interaction Di Felice42 force others (e.g., Koch and Hill43)

model A

model B

refs 3, 6-13 refs 14-21 refs 22-30 refs 31-38 refs 3, 6-8, 12-21, refs 33, 37, 38 24, 26, 28-30 ref 25 refs 4, 5, 31, 32, 34-36 refs 9, 22, 23, 27, 39

study various gas-solid flow phenomena, particularly in fluidization as listed in Table 1. However, as also listed in Table 1, differences exist in the implementation of the CCDM. The differences can be categorized according to three aspects: governing equations, schemes of coupling the gas and solid phases, and correlations to calculate the particle-fluid interaction force. Two sets of governing equations, referred to as models A and B1, are used in the TFM. Correspondingly, there are also two sets of governing equations or two model formulations in the CCDM,44 although they describe the solid phase at an individual-particle level instead of at a localaverage, computational-cell level. Different schemes have been used to couple the gas and solid phases through a particle-fluid interaction force, namely, from particle scale to computational-cell scale or the opposite. In addition, different correlations have been used to calculate the particle-fluid interaction force. Although these aspects are probably difficult to fully resolve because they might be related to the limit of the current theoretical and experimental understanding, it is important at least to know their effects on simulation results.

10.1021/ie049387v CCC: $27.50 © 2004 American Chemical Society Published on Web 11/18/2004

Ind. Eng. Chem. Res., Vol. 43, No. 26, 2004 8379

This paper presents a study of differences in the above aspects. The origins of problem are discussed first. Then, emphasis is placed on the effects of the model formulations. For the purpose of illustration, gas fluidization of particle mixtures is simulated with both model formulations, facilitated by different correlations to calculate the particle-fluid interaction force. The results are compared in terms of flow pattern and mixing/ segregation kinetics and explained in terms of the particle-fluid interaction force. Physical experiments are also conducted under comparable conditions to generate a sound basis for assessing the treatments associated with the model formulations. Theoretical Consideration Governing Equations. In the CCDM, the gas phase is treated as a continuous phase, and its flow complies with the law of conservation of mass and momentum in terms of local-average variables. The governing equations for the gas phase are actually the same as those used in the TFM, for which two formulations, referred to as models A and B, have been proposed.1 Model A assumes that the pressure drop is shared between the gas and solid phases, and model B assumes that it applies to the gas phase only. Both formulations have been used in the CCDM, and their corresponding governing equations for the gas phase can be expressed as follows

Conservation of mass ∂g + ∇‚(gug) ) 0 ∂t

(1)

Conservation of momentum model A ∂(Fggug) + ∇‚(Fggugug) ) ∂t -g∇P - FA + ∇‚(gτ) + Fggg (2a) model B ∂(Fggug) + ∇‚(Fggugug) ) ∂t -∇P - FB + ∇‚(gτ) + Fggg (2b) where u and P are the fluid velocity and pressure, respectively, and τ, g, Fg, and ∆V are the fluid viscous stress tensor, fluid volume fraction, fluid density, and computational cell volume, respectively. FA and FB are the volumetric particle-fluid interaction forces for the two models. Note that the definition of the total particlefluid interaction force can vary depending on how one interprets eq 2. Nonetheless, it is well-established that there is a link between FA and FB, given by FB ) FA/g - Fgsg.1 The solid phase is treated by a discrete approach, involving two types of motion: translational and rotational. The forces acting on a particle include gravity, contact forces among particles and between particles and the wall, and interaction forces between particles and the fluid. Other forces such as the van der Waales force associated with fine particles and the capillary force associated with wet particles can also be included.5,17-19,35,36 For simplicity, they are not considered

in this work. Therefore, the governing equations for particle i can be expressed as follows

Translational motion dup,i

mi

dt

ki

) fp-f,i +

(fc,ij + fd,ij) + Fp,iVp,ig ∑ j)1

(3a)

Rotational motion Ii

dωi dt

ki

)

Tij ∑ j)1

(3b)

where mi, Ii, ki, up,i, and ωi are the mass, moment of inertia, number of contacting particles, and translational and rotational velocities of particle i, respectively, and fc,ij, fd,ij, and Ti,j are the contact force, viscous contact damping force, and torque, respectively, between particles i and j. These interparticle forces and torques are summed over the ki particles in contact with particle i. fp-f,i is the total particle-fluid interaction force on a particle, which includes many forms, such as the fluid drag force fdrag,i, the buoyancy force, the lifting force, the virtual mass force, the Basset force, etc.45 For simplicity, consistent with the traditional TFM treatment of gas fluidization, only the buoyancy and fluid drag forces are considered in this work. The particlefluid interaction force can be determined from the empirical correlations discussed later. Corresponding to eq 2 for the fluid phase, two approaches can be used to determine the particle-fluid interaction force. For convenience, they are also referred to as models A and B, although they are obtained at different time and length scales.44,46 Thus, we have

Model A

fp-f,i ) -Vp,i∇Pi + fA

(4a)

Model B

fp-f,i ) FgVp,ig + fB

(4b)

As for the gas phase, there is a relation between fA and fB, given by fB ) fA/g. Models A and B for the gas and solid phases must correctly match in order to generate correct quantitative results, as noted by Xu and Yu.44 Different forms of the particle-fluid interaction force can be identified between models A and B. In model A, one part is the buoyancy force related to the total pressure drop ∇Pi, and the other part, fA, is the fluid drag force multiplied by fluid volume fraction gfdrag,i. In contrast, in model B, the particle-fluid interaction force comprises the buoyancy force, which is related only to the static pressure drop ∇P0, and the fluid drag force fdrag,i. The total pressure drop ∇P mainly consists of three parts: the hydrostatic pressure drop ∇P0 due to the gravity of the gas, the dynamic pressure drop or the so-called piezometric or manometric pressure drop ∇Pd due to the relative motion between the gas and the particles, and the pressure drop due to friction between the gas and wall ∇Pw. ∇Pw is often ignored. Thus, eqs 4a and 4b can be rewritten as

Model A

fp-f,i ) -Vp,iFgg + Vp,i∇Pd,i + gfdrag,i (5a)

Model B

fp-f,i ) -Vp,iFgg + sfdrag,i + gfdrag,i

(5b)

For a system composed of monosized particles that is uniform and has no acceleration in either phase, the

8380

Ind. Eng. Chem. Res., Vol. 43, No. 26, 2004

pressure drop is just the static pressure drop ∇P0 and the dynamic pressure drop ∇Pd. The pointwise values around every particle, ∇ Pd,i, are identical and can be replaced with the local-average value ∇Pd. The particle volume Vp,i can be expressed as s/n, where n is the number of particles per unit volume and the relation between the dynamic pressure drop and drag force acting on a particle is ∇Pd ) nfdrag,i. Therefore, the second term on the right-hand side of eq 5a, Vp,i∇Pd,i, can be expressed as (s/n)nfdrag,i ) sfdrag,i. Thus, eqs 5a and 5b produce the same results for the solid phase. That is, although different in form, models A and B are the same for this simple system. However, this consideration is not applicable to fluidization for the following reasons: The pressure drop used to calculate the particle-fluid interaction force in model A should be based on the pointwise value around each particle i, ∇Pd,i. However, the only pressure drop one can obtain from a CCDM simulation is the localaverage value ∇Pd as a result of the nature of the continuum approach. In fluidization, particles are not uniformly distributed in space and have different velocities in both magnitude and direction, as is the case for the gas phase. Consequently, the pointwise value ∇Pd,i does not equal the local-average value ∇ Pd. This difference results in differences between the models A and B simulations, as demonstrated in this work. Coupling between the Discrete Solid Phase and the Continuum Gas Phase. In gas fluidization, there are strong interactions between the continuum gas and the discrete particles. The coupling between the two phases is achieved through the particle-fluid interaction force, which is at the local-average or computationalcell level for the gas phase (FA or FB) and at the individual-particle level for the solid phase (fAor fB). Different schemes have been proposed to couple the two phases modeled at different length scales. According to scheme 1, the force from the particles to the gas phase is calculated by a local-average method as used in the TFM, whereas the force from the gas phase to each particle is calculated separately according to individual-particle velocity. According to scheme 2, the force from the particles to the gas phase is calculated first at a local-average scale as used in scheme 1. This value is then distributed to individual particles according to a certain average rule. Typically, for monosized particles, the average rule used is given by

f)

F∆V kc

(6)

where kc is the number of particles in a computational cell, F is the volumetric particle-fluid interaction force, ∆V is the volume of the computational cell, and f is the force at the particle scale. F and f can be in the form of model A or model B. According to scheme 3, at each time step, the particlefluid interaction forces on individual particles in a computational cell are calculated first, and the values are then summed to produce the particle-fluid interaction force at the cell scale, typically given by kc

F)

f ∑ i)1 ∆V

(7)

According to the Newton’s third law of motion, the force of the solid phase acting on the gas phase should be equal to the force of the gas phase acting on the solid phase but in the opposite direction. Scheme 1 does not guarantee that this condition can always be satisfied. Consequently, it is not reasonable. Indeed, this scheme was used only in the early stage of CCDM development,3,7 although it can still be found occasionally in the literature as listed in Table 1. Scheme 2 can satisfy Newton’s third law. However, it uniformly distributes the interaction force among the particles in a computational cell irrespective of the different behaviors of these particles in the cell. This scheme cannot fully represent reality, as the particle-fluid interaction forces for the particles in the cell should differ for nonuniform gas-solid flow. In addition, in the calculation of the particle-fluid interaction force F, a mean particle velocity has to be used. The appropriate method for calculating this mean particle velocity is still an open question, particularly for multisized particle systems. Scheme 3 can overcome the above problems associated with schemes 1 and 2. Indeed, this scheme has been widely accepted since its first introduction by Xu and Yu,4 as seen in Table 1. Particle-Fluid Interaction Force. The particlefluid interaction force, mainly in the form of a fluid drag force for gas-solid flow systems, is often the driving force in gas fluidization. The accurate representation of the fluid drag force is important for a quantitative simulation and analysis of a particulate system. Because of the complexity of the problem, a precise prediction of this force is not available. Only for an isolated particle under low Reynolds numbers (Stokes’ region) is an analytical solution of the Navier-Stokes equations possible.47 To date, understanding the effects of different arrangements of surrounding particles on the drag coefficient for a considered particle is still an active research area (see ref 48, for example). Recently, numerical techniques have been developed to quantify the drag force through a detailed study of the forces from a fluid on a particle, such as direct numerical simulation (DNS)49 and the latticeBoltzmann (LB) method.43,50 The combination of such a method and the DEM, e.g., DNS-DEM and LBDEM, would offer a rational approach not only to determining the particle-fluid interaction forces but also to modeling of the particle-fluid flow without the problems considered in this work.51 However, limited by the current computational capability, to date, such studies have been applied only to very simple systems.52-54 Therefore, at this stage of development, the calculation of fA or fB is made mainly based on the existing correlations for pressure drop, ∇P, in the literature. This consideration also explains why we still use the concepts of models A and B, unlike Kafui et al.25,55 The modeling of fluid flow in the CCDM is still at a local-average, computational-cell level, built on the basis of the TFM, rather than at the subparticle or even particle scale.5,51 For a system uniform for both the gas and solid phases, the pressure drop is uniformly shared by the particles. For such a system, the following relations can be derived: for model B

fB )

∇Pd π 3 1 ) di ∇Pd n 6 s

(8)


and for model A

Table 2. Parameters Used in the Simulations

g π g f ) ∇Pd ) di3 ∇Pd n 6 s A

parameter

(9)

Packed beds of monosized particles with uniform fluid flow are typical uniform systems, as discussed earlier. Therefore, the above equation offers a simple and convenient way to determine the particle-fluid interaction force on a particle scale, as used by Di Felice.42 Many correlations have been proposed in the literature to determine the pressure drop of such a packed bed. In principle, all can be used for this purpose. In the CCDM application, as seen in Table 1, probably because of the inheritance from the previous TFM work,1 the most common form is the combination of the Ergun equation40 for fluid volume fractions less than 0.8 and Wen and Yu’s correlation41 for fluid volume fractions greater than 0.8. This treatment causes a discontinuous drag force at a fluid volume fraction of 0.8. Hence, the Di Felice correlation42 is becoming increasingly popular. In addition to these two forms, other forms43,56,57 have also been used. To date, all equations in the CCDM simulation have been derived for monosized particle systems. Strictly speaking, the appropriate method of describing the fluid drag force for multisized particle systems is still an open research topic. At present, the correlations for monosized particle systems are directly applied to multisized particle systems.22,24,32 Note that Hoomans et al.24 and Bokkers et al.22 used a hard-sphere model in their DPM, which focused on the kinematic behavior of the particles and did not explicitly consider the interaction forces between particles. On the other hand, Feng et al.32 used a three-dimensional DEM or soft-sphere model for the solid phase. Their approach is employed in this work, as described below. Results and Discussion The above discussion clearly indicates that there are three important issues in implementing the CCDM approach to simulate gas-solid flow. The issue of the coupling scheme has been well clarified in the literature, and scheme 3 should be adopted. The issue of calculating the particle-fluid interaction force is actually related to the accuracy of the correlations proposed in the literature. As noted by Feng and Yu,46 it might be difficult to do much about this matter. The fact here is that there is a large error associated with such a correlation. For example, the Ergun equation, the most popular approach for calculating the pressure of a bed, can be expected to predict experimental results with an accuracy of only (50%, as discussed by Macdonald et al.58 The CCDM simulations with different correlations, particularly for those relatively well received, might simply produce results with differences largely within this error range. Indeed, it has been reported that the Di Felice correlation and the combined Ergun and Wen and Yu correlation yield similar CCDM results.27,46 It should be noted that, within the framework of the TFM, the effects of these and other correlations on the computed fluidization behavior have been examined by various investigators.59-61 Therefore, the issue of model formulation is the most critical in CCDM development, which is the focus of this work. For monosized particle systems, its effects, in terms of models A and B, on the simulated behavior

units

Particle Phase density kg m-3 Young’s modulus N m-2 Poisson ratio sliding friction coefficient damping coefficient particle diameter m flotsam jetsam number of particles flotsam jetsam time step s viscosity density CFD cell width height bed geometry width height thickness gas distributor time step

Gas Phase kg m-1 s-1 kg m-3

value 2500 1.0 × 108 0.3 0.3 0.2 0.001 0.002 22 223 2777 2.5 × 10-6 1.8 × 10-5 2500

m m

0.00325 0.00325

m m m

0.065 0.26 0.0081 uniform 2.5 × 10-6

s

have been studied by the TFM and the CCDM. Bouillard et al.59 reported that there is no significant difference between models A and B simulations in the TFM approach. Feng and Yu46 recently obtained a similar conclusion for the CCDM approach. However, this consideration does not apply to the fluidization of particle mixtures. In this work, the gas fluidization of binary mixtures of particles is simulated to highlight this point. Simulation Conditions. The fluidization of a binary mixture (flotsam diameter ) 0.001 m and jetsam diameter ) 0.002 m) in a rectangular bed is considered in this work. In total, 25 000 particles are used in a simulation with particles of each size accounting for 50% by weight. Table 2 lists the parameters used in this work. The equations used to calculate the interaction forces among particles and between particle and fluid were provided in our previous work32 and are therefore not repeated here. Unless otherwise stated, in this work, the Di Felice correlation is used for the calculation of the particle-fluid interaction force. The correlation used in the work of Feng et al.32 is in the form of model B. The form for model A can be obtained according to the relationship fA ) gfB as discussed earlier. In our previous study,32 periodic boundary conditions were applied to the front and rear walls to present the full 3D motion of the particles with a relatively small number of particles. However, this causes a difficulty for direct comparisons with experimental results because the front and rear walls are always needed to support particles in physical experiments. Therefore, in this work, the front and rear walls are included in all simulations. The physical properties of walls are assumed to be the same as those of the particles. The flow of the gas is still assumed to be two-dimensional, achieved by using a bed with a relatively small thickness. Each simulation is started with the random generation of particles according to their preset proportions without overlaps in the bed, followed by a gravitational settling process for 0.6 s. This produces a well-mixed

8382


Figure 1. Particle flow patterns at different times for model A when the superficial gas velocity is (a) 1.2, (b) 1.3, (c) 1.4, and (d) 1.6 m/s.

packed bed that is used as the initial state for all fluidization cases considered. Gas is then injected uniformly at the bottom of the bed to fluidize the bed at a preset velocity. The gas velocities used in this study are 1.2, 1.3, 1.4, and 1.6 m/s. The conditions for all computations are exactly the same except for the aspect to be examined. These include the same initial packing, material properties, CCDM formulations, and numerical algorithms. This treatment is important for generating comparable results.32 Model A vs Model B Simulations. The simulated results were first checked in terms of flow patterns. For better visualization, only particles whose center points are between 1.5 and 2.5 mm in the thickness direction are shown. Figure 1 shows the solid flow patterns for the model A simulation at different velocities. When the gas injection velocity is 1.2 m/s (Figure 1a), once gas is introduced into the bed, size segregation gradually

Figure 2. Particle flow patterns at different times for model B when the superficial gas velocity is (a) 1.2, (b) 1.3, (c) 1.4, and (d) 1.6 m/s.

appears, with flotsam moving upward and remaining in the fluidized state at the top layer. The decrease of the concentration of flotsam particles in the bottom part causes an increase of the minimum fluidization velocity of the particles there and gradually leads to the formation of a defluidized layer. In this layer, some flotsam particles are entraped by surrounding jetsam particles and cannot move upward. Finally, a macroscopically stable state is reached. Two layers can be clearly identified: the top layer, rich in flotsam and in a fluidized state, and the bottom layer, rich in jetsam and in a defluidized state. With increasing gas velocity, the bed becomes more fluidized, giving different segregation patterns as shown in parts b-d of Figure 1. When the gas velocity is 1.6 m/s, the bed is almost fully fluidized, producing a resonably well-mixed state, with a small number of jetsam particles sitting on the bottom.


Figure 3. Variation of pressure drop with time for different models (left, model A; right, model B) when the superficial gas velocity is (a) 1.2 and (b) 1.6 m/s.

The same description applies to the model B simulations as shown in Figure 2. However, a comparison between the models A and B simulations indicates that model B produces a more significant segregation and a clearer interface between the fluidized and defluidized layers at a given velocity. The difference is less significant at high gas velocities, e.g., 1.6 m/s, when the whole bed is fluidized. As with the solid flow pattern, there is a difference in the gas pressure drop. Figure 3 shows the bed pressure drop measured at the bottom of the bed for two gas velocities. It is obvious that, when the gas velocity is 1.2 m/s, model A produces a larger fluctuation than model B. The model B simulation results in a larger defluidized layer than the model A simulation and, hence, less fluctuation in the pressure drop. The fluctuation in the pressure drop increases with increasing gas velocity, but the difference between the two models becomes less obvious once the bed is fully fluidized. Comparison is further made by quantifying the mixing or segregation kinetics in terms of a mixing index. The Lacey mixing index62 is used for this study. The method of calculation for both models is the same as that used in our previous work.32 Figure 4 shows the variation of the mixing index with time at different gas velocities. Starting at the well-mixed packing state, the mixing index decreases as size segregation develops. In general, at a given velocity, the decreasing rate for model B is higher than that for model A, and the index reaches a lower value for model B when the macroscopically stable state is reached. Physical vs Numerical Experiments. The differences between models A and B described in the above section offer an excellent opportunity to assess the relative validity of the two formulations. For this purpose, physical experiments were conducted under conditions comparable to those used in the simulations. The walls of the experimental cell were made of Perspex sheet. Glass beads of diameters 1 and 2 mm were used as the packing/fluidized materials. The rest of the parameters, including the bed geometry and gas flow conditions, were exactly the same as used in the

Figure 4. Mixing index as a function of time for different model formulations (black, model A; gray, model B) when the gas injection velocity is (a) 1.2, (b) 1.3, (c) 1.4, and (d) 1.6 m/s.

numerical simulations. The process was recorded by a digital video camera. Different sized particles can be identified by the colors reflected through light (dark for flotsam and bright for jetsam). Figure 5 shows the solid flow patterns measured. A comparison between the physical and numerical experiments (Figures 1, 2, and 5) indicates that both model formulations are capable of producing reasonable flow patterns, at least qualitatively, but generally speaking, the model B results are more consistent with the experimental observations. It is important to quantify the difference between the measured and simulated results. However, such quantitative analysis is difficult to achieve with the present experimental techniques. Thus, the present quantitative analysis is made for “frozen” beds only, obtained by stopping the gas flow and allowing the particles to fall under the force of gravity to form a packed bed. Segregation patterns in the frozen beds were then quantified by analyzing the concentration of the layered samples taken from the top to the bottom of the bed. The same operation was applied to the numerical simulation as well.

8384


Figure 5. Solid flow patterns observed in the experiment when the gas injection velocity is (a) 1.2, (b) 1.3, (c) 1.4, and (d) 1.6 m/s.

The experimental results might be affected by operation. To check this effect, the same experiment was conducted three times. Figure 6 shows the concentration profile of flotsam along the vertical direction for three runs when the gas injection velocity was 1.3 m/s and after the macroscopically stable state had been reached. Clearly, the experimental results are quite reproducible. Repetitive tests were not conducted for the numerical experiments because they are always reproducible un-

less there are changes in simulation conditions. However, the profile might vary with time because of the transient behavior of fluidization. This aspect was examined by analyzing the results obtained when the gas supply was shut off at different times, as shown in Figure 7. The results indicate that, even though bubbling can cause transient solid flow pattern, such bubbling does not result in much difference in the segregation pattern of the frozen bed.


Figure 6. Concentration of flotsam along the vertical direction, measured from the physical experiment for three runs when the gas injection velocity is 1.3 m/s.

Figure 7. Concentration of flotsam along the vertical direction, determined from the numerical experiment with models A and B when the bed is frozen at different times and the gas injection velocity is 1.3 m/s.

Figure 8 shows the measured and simulated profiles along the vertical direction at different velocities, all obtained when the beds were in the macroscopically stable state. Generally, two layers can be clearly seen at low velocities (v ) 1.2, 1.3, and 1.4 m/s): a bottom layer rich in jetsam and a top layer rich in flotsam. At these low velocities, the model B simulation shows better agreement with the physical experiment than the model A simulation (Figure 8a-c). At higher velocity (v ) 1.6 m/s), when the bed is fully fluidized, simulations with both models A and B agree with the experiments reasonably well. In fact, referring to the results in Figures 6 and 7, the difference between the model A

simulation and the experimental observations is almost within the experimental error. Role of the Particle-Fluid Interaction Force. The flow and segregation of the particles in a fluidized bed is controlled by the interactions among particles, between particle and fluid, and between particle and wall, in addition to the gravitational force. An analysis of these interaction forces and their spatial and temporal evolution can provide a better understanding of the underlying mechanisms as demonstrated in our previous work.5,32 Here, we extend that analysis to investigate the effects of different model formulations. Because the same model for particle-particle interactions is employed, the difference in flow behavior and mixing kinetics between models A and B must be mainly associated with the different particle-fluid interaction forces from the two model formulations. Therefore, the present analysis is focused on the particle-fluid interaction force. The particle-fluid interaction forces for models A and B are calculated according to eqs 5a and b, respectively. In model A, the pressure gradient should be based on the pointwise value around each particle i. However, this pointwise value is not available, and hence, it is replaced by the local-average pressure drop. This treatment was widely used in the previous model A simulations as listed in Table 1. In model B, the pressure drop relates to only the hydrostatic pressure drop; it does not require this approximation. To examine their effects on fluid drag force, we first applied the two model formulations to a well-mixed packed bed composed of the same particle mixture, i.e., the initial packing for fluidization in this work. Figure 9 shows the particle-fluid interaction forces calculated with the two model formulations as a function of the fluid volume fraction for a superficial gas velocity of 1.0 m/s. For better comparability and representativeness, the interaction force on a particle is expressed as a dimensionless force by dividing the particle weight. Thus, as shown in Figure 9, a 1-mm particle has a larger value than a 2-mm particle in both model A and model B, which explains why the small particles (flotsam) are driven to the top layer once gas injection starts. More importantly, it shows that model B gives a larger

Figure 8. Physical vs numerical experiments in terms of the concentration of flotsam along the vertical direction when the gas injection velocity is (a) 1.2, (b) 1.3, (c) 1.4, and (d) 1.6 m/s.

8386


Figure 9. Particle-fluid interaction force at different porosities when the superficial velocity is 1.0 m/s: line 1, model A with the dynamic pressure shared according to volume fraction, d ) 2 mm; line 2, model A with the dynamic pressure shared according to volume fraction, d ) 1 mm; line 3, model B, d ) 2 mm; line 4, model B, d ) 1 mm; line 5, model A with the dynamic pressure shared according to surface fraction, d ) 2 mm; line 6, model A with the dynamic pressure shared according to surface fraction, d ) 1 mm.

difference between 1- and 2-mm particles than model A. This must be a major factor responsible for the different flow patterns between the models A and B simulations as observed in Figures 1-4. To confirm this conclusion, the particle-fluid interaction force in the numerical simulation was further checked. Because the segregation between flotsam and jetsam particles mainly happens in the vertical direction, the present analysis focuses only on the forces in this direction and, for simplicity, their respective mean values. Figure 10 shows the variation of the mean particle-fluid interaction force with time for a gas injection velocity of 1.3 m/s. Again, the force acting on a particle is made dimensionless by dividing the particle weight. The upward particle-fluid interaction force is larger than the weight for flotsam and less than the weight for jetsam. This provides a driving force for segregation. As discussed by Feng et al.,32 two stages can be identified: the transient stage and the stable stage. In the transient stage, the value of the mean particle-fluid interaction force decreases for flotsam and increases for jetsam. The difference in the particle-fluid interaction force between flotsam and jetsam is larger for model B than for model A, which agrees well with the theoretical considerations, as shown in Figure 9. In the stable stage, the values simply fluctuate around their respective means. The difference observed between Figures 9 and 10 is mainly associated with the approximation involved in model A. The total pressure drop here includes the static pressure drop, ∇P0, due to the gravity of the gas and the dynamic pressure drop, ∇ Pd, due to the particlefluid interaction. The buoyancy force related to the static pressure drop ∇P0 is distributed to individual particles according to their volumes. For the dynamic pressure drop, though, it should be based on the pointwise value ∇ Pd,i, which is not available and replaced by localaverage value in a simulation. Thus, a difference between models A and B results. Modification of the Model A Approach. There are two methods to remedy the inadequacy in distributing the dynamic pressure to individual particles. The first one is to use the pointwise value ∇Pd,i rather than the local-average value ∇Pd. Obviously, this pointwise value is impossible to obtain within the present CCDM model

Figure 10. Variation of the mean particle-fluid interaction forces on the flotsam and jetsam with time in the vertical direction when the gas injection velocity is 1.3 m/s: (a) model A, (b) model B.

Figure 11. Particle flow patterns at different times when the gas injection velocity is 1.3 m/s for model A with the dynamic pressure distributed to particles according to their surface fraction.

framework, but it can be treated as the sum of contributions from pointwise particles and then redistributed to the particles. In this case, it equals sfdrag,i. Thus, this treatment is actually identical to the treatment applied in model B. The second method is to continue to use the local-average pressure drop ∇Pd, but to calculate the hydrodynamic force acting on the particles in a cell as s∇Pd and to assign this value to each particle according to its surface area rather than volume. This treatment, although still an approximation, should be more reasonable, because the dynamic pressure drop comes from the particle-fluid interaction on the surface of a particle. The second method was tested in this work. As shown in Figure 9, the particle-fluid interaction force for model A based on surface area produces a similar value to the model B calculation. The test was then extended to the fluidization of the particle mixture considered.


Figure 13. Mixing index as a function of time for different models when the superficial gas velocity is 1.3 m/s: lines 1 and 4, model A with the dynamic pressure shared according to volume fraction; lines 2 and 5, model A with the dynamic pressure shared according to surface fraction; lines 3 and 6, model B. The Ergun and Wen and Yu correlations are used for calculating the fluid drag force for lines 1-3, and the Di Felice correlation is used for lines 4-6.

Figure 12. Particle flow patterns at different times when the gas injection velocity is 1.3 m/s for different models in which the fluid drag force is calculated according to a combination of the Ergun and Wen and Yu correlations: (a) model A with the dynamic pressure shared according to volume fraction, (b) model A with the dynamic pressure shared according to surface fraction, (c) model B.

The gas injection velocity used was 1.3 m/s. At each time step, the hydrodynamic force in a computational cell was distributed to the particles according to their fractional contributions to the total surface area in the cell. The solid flow patterns at different times are shown in Figure 11, which indeed match reasonably well the model B results in Figure 2b. The above results highlight the need to distinguish between the dynamic pressure and the static pressure. The former results from the interactions between particles and fluid on the surface of the particles and, hence, should be distributed to particles according to their contributions to the surface area in a computational cell. The latter comes from the volume replacement of the particles in a fluid, and hence, its distribution to particles should be based on their contributions to the replaced volume. This point is very important as it could represent the origin leading to the differences between models A and B. For monosized particles, there is no difference in surface and volume fractions. Consequently, for such particle systems, there is no significant difference between the models A and B simulations.46 However, minor differences can always be observed. Such differ-

Figure 14. Concentration of flotsam along the vertical direction when the gas velocity is 1.3 m/s: lines 1 and 4, model A with the dynamic pressure shared according to volume fraction; lines 2 and 5, model A with the dynamic pressure shared according to surface fraction; lines 3 and 6, model B. The Ergun and Wen and Yu correlations are used for calculating the fluid drag force for lines 1-3, and the Di Felice correlation is used for lines 4-6.

ences stem from the difference in the scheme for coupling gas and solid phases. In the model A approach, the pressure drop at a computational-cell level has to be distributed to individual particles. The treatment depends not only on the surface areas of the particles but also on the relative velocity between the particles and the fluid. Often, the latter is simply not considered, as is the case in all of the previously reported CCDM simulations. In other words, the difference between models A and B is probably inevitable. Information transfer from a large scale to a small scale is always problematic. From this point of view, model B is favored as it does not involve as many such data transfers as model A. One might question whether the above remarks would be incorrect if different correlations had been used to calculate the particle-fluid interaction force. To answer this question, simulations were conducted for models A and B, facilitated with the Ergun and Wen and Yu correlations. Figure 12 shows the flow patterns of different models when the gas injection velocity is 1.3 m/s. It can be seen that the difference does exist for the models A and B simulations (part a vs part c of Figure

8388


12). When model A is modified by sharing the dynamic pressure drop according to the surface fraction, it produces flow patterns in better agreement with those from the model B simulation (part b vs part c of Figure 12). This can be further quantitatively confirmed by examination of the mixing kinetics and the concentration profiles determined after the macroscopically steady state has been reached, as shown in Figures 13 and 14. A comparison of the results in the two figures also indicates that different correlations for calculating the particle-fluid interaction force can produce different results. However, the difference is not as significant as that from different model formulations. Conclusions The different implementations of the CCDM in the literature have been summarized and discussed in terms of the governing equations, the coupling schemes between the continuum gas and discrete solid phases, and the correlation for calculating the particle-fluid interaction force. Numerical experiments for the fluidization of binary mixtures of particles show that there is a significant difference between the models A and B simulations, reflected in the solid flow pattern and mixing/segregation kinetics. The difference is shown to stem from the different particle-fluid interaction forces of the two models. It is important to distinguish the static pressure drop and dynamic pressure drop, which are related to particle volume and surface area, respectively. The difference between the models A and B formulations mainly results from the distribution of the force associated with the dynamic pressure drop at a local-average scale to individual particles in model A. It can be reduced by modifying the method of distribution according to the surface area rather than the volume of the particles. Different correlations for calculating the particlefluid interaction force can produce simulated results that are quantitatively different but qualitatively similar. Considering the large error associated with the development of these models, such differences are probably not as critical as the model formulation in producing consistent results. Comparison with the physical experiments conducted under comparable conditions and analysis of the numerical scheme for implementing a CCDM simulation suggest that the model B formulation is favorable. Further studies are necessary to fully clarify this important issue. Acknowledgment The authors are grateful to the Australia Research Council (ARC) for financial support, the Australian Partnership for Advanced Computing (APAC) for computing support, and the Australian Centre for Advanced Computing and Communications (AC3) for the use of their computational facilities. DR. B. H. Xu of the University of Leeds (U.K.) is also thanked for the beneficial discussion at the early stage of this work. Finally, Y.Q.F. is grateful to the University of New South Wales for providing a scholarship. Nomenclature d ) particle diameter, m f ) force, N F ) volumetric fluid-particle interaction force, N m-3

g ) gravitational acceleration, m s-2 I ) rotational inertial momentum of a particle, kg m2 kc ) number of particles in a computational cell ki ) number of particles in contact with particle i m ) mass, kg n ) number of particles per unit volume P ) pressure, Pa ∇P0 ) hydrostatic pressure drop, Pa m-1 ∇Pd ) hydrodynamic pressure drop, Pa m-1 t ) time, s T ) torque, N m u ) velocity, m s-1 V ) volume, m3 ∆V ) volume of a computational cell, m3 Greek Letters ) volume fraction µ ) viscosity, kg m-1 s-1 F ) density, kg m-3 τ ) continuum-phase viscous stress tensor, kg m-1 s-2 ω ) rotational velocity of a particle, s-1 Subscripts c ) contact d ) damping f ) fluid phase g ) gas i ) particle i ij ) between particles i and j p ) particle p-f ) between particle and fluid s ) solid phase Superscripts A ) model A B ) model B

Literature Cited (1) Gidaspow, D. Multiphase Flow and Fluidization; Academic Press: San Diego, 1994. (2) Cundall, P. A.; Strack, O. D. L. A discrete numerical model for granular assemblies. Geotechnique 1979, 29, 47. (3) Tsuji, Y.; Kawaguchi, T.; Tanaka, T. Discrete particle simulation of two-dimensional fluidised bed. Powder Technol. 1993, 77, 79. (4) Xu, B. H.; Yu, A. B. Numerical simulation of the gas-solid flow in a fluidised bed by combining discrete particle method with computational fluid dynamics. Chem. Eng. Sci. 1997, 52, 2786. (5) Yu, A. B.; Xu, B. H. Particle scale modelling of particlefluid flow in fluidization, J. Chem. Technol. Biotechnol. 2003, 78, 111. (6) Bin, Y.; Zhang, M. C.; Dou, B. L.; Song, Y. B.; Wu, J. Discrete particle simulation and visualized research of the gas-solid flow in an internally circulating fluidized bed. Ind. Eng. Chem. Res. 2003, 42, 214. (7) Hoomans, B. P. B.; Kuipers, J. A. M.; Briels, W. J.; van Swaaij, W. P. M. Discrete particle simulation of bubble and slug formation in a two-dimensional gas-fluidised bed: A hard-sphere approach. Chem. Eng. Sci. 1996, 51, 99. (8) Limtrakul, S.; Chalermwattanatai, A.; Unggurawirote, K.; Tsuji, Y.; Kawaguchi, T.; Tanthapanichakoon, W. Discrete particle simulation of solids motion in a gas-solid fluidized bed. Chem. Eng. Sci. 2003, 58, 915. (9) Yuu, S.; Umekage, T.; Johno, Y. Numerical simulation of air and particle motions in bubbling fluidized bed of small particles. Powder Technol. 2000, 110, 158. (10) Kawaguchi, T.; Tsuji, Y. Numerical simulation of twodimensional fluidized beds using the discrete element mehod (comparison between the two-and three-dimensional models). Powder Technol. 1998, 96, 129.

Ind. Eng. Chem. Res., Vol. 43, No. 26, 2004 8389 (11) Kawaguchi, T.; Sakamoto, M.; Tanaka, T.; Tsuji, Y. Quasithree-dimensional numerical simulation of spouted beds in a cylinder. Powder Technol. 2000, 109, 3. (12) Ouyang, J.; Li, J. Particle-motion-resolved discrete model for simulating gas-solid fluidization. Chem. Eng. Sci. 1999, 54, 2077. (13) Ouyang, J.; Li, J. Discrete simulation of heterogeneous structure and dynamic behavior in gas-solid fluidization. Chem. Eng. Sci. 1999, 54, 5427. (14) Han, T.; Levy, A.; Kalman, H. DEM simulation for attrition of salt during dilute-phase pneumatic conveying. Powder Technol. 2003, 129, 92. (15) Kaneko, Y.; Shiojima, T.; Horio, M. DEM simulation of fluidized beds for gas-phase olefin polymerization. Chem. Eng. Sci. 1999, 54, 5809. (16) Kuwagi, K.; Mikami, T.; Horio, M. Numerical simulation of metallic solid bridging particles in a fluidised bed at high temperature. Powder Technol. 2000, 109, 27. (17) Mikami, T.; Kamiya, H.; Horio, M. Numerical simulation of cohesive powder behavior in a fluidized bed. Chem. Eng. Sci. 1998, 53, 1927. (18) Rhodes, M. J.; Wang, X. S.; Ngyen, M.; Stewart, P.; Liffman, K. Use of discrete element method simulation in studying fluidization characteristics: Influence of interparticle force. Chem. Eng. Sci. 2001, 56, 69. (19) Rhodes, M. J.; Wang, X. S.; Nguyen, M.; Stewart, P.; Liffman, K. Study of mixing in gas-fluidized beds using a DEM model. Chem. Eng. Sci. 2001, 56, 2859. (20) Rong, D.; Mikami, T.; Horio, M. Particle and bubble movements around tubes immersed in fluidized bedssA numerical study. Chem. Eng. Sci. 1999, 54, 5737. (21) Wang, X. S.; Rhodes, M. J. Determination of particle residence time at the walls of gas fluidized beds by discrete element method simulation. Chem. Eng. Sci. 2003, 58, 387. (22) Bokkers, G. A.; van Sint Annaland, M.; Kuipers, J. A. M. Mixing and segregation in a bi-disperse gas-solid fluidized bed: A numerical and experimental study. Powder Technol. 2004, 140, 176. (23) Helland, E.; Occelli, R.; Tadrist, L. Numerical study of cluster formation in a gas-particle circulating fluidized bed. Powder Technol. 2000, 110, 210. (24) Hooman, B. P. B.; Kuipers, J. A. M.; van Swaaij, W. P. M. Granular dynamics simulation of segregation phenomena in bubbling gas-fluidised beds. Powder Technol. 2000, 109, 41. (25) Kafui, K. D.; Thornton, C.; Admas, M. J. Discrete particlecontinuum fluid modelling for gas-solid fluidised beds. Chem. Eng. Sci. 2002, 57, 2395. (26) Li, J.; Kuipers, J. A. M. Effect of pressure on gas-solid flow behavior in dense gas-fluidized beds: A discrete particle simulation study. Powder Technol. 2002, 127, 173. (27) Li, J.; Kuipers, J. A. M. Gas-particle interactions in dense gas-fluidized beds. Chem. Eng. Sci. 2003, 58, 711. (28) van Wachem, B. G. M.; van der Schaaf, J.; Schouten, J. C.; Krishna, R.; van den Bleek, C. M. Experimental validation of Lagrangian-Eulerian simulation of fluidized beds. Powder Technol. 2001, 116, 155. (29) van Wachem, B. G. M.; Almstedt, A. E. Methods for multiphase computational fluid dynamics. Chem. Eng. J. 2003, 96, 81. (30) Zhou, H.; Flamant, G.; Gauthier, D.; Lu, J. Lagrangian approach for simulating the gas-particle flow structure in a circulating fluidized bed riser. Inter. J. Multiphase Flow 2002, 28, 1801. (31) Feng, Y. Q.; Pinson, D.; Yu, A. B.; Chew, S. J.; Zulli, P. Numerical study of gas-solid flow in the raceway of a blast furnace. Steel Res. Int. 2003, 74, 523. (32) Feng, Y. Q.; Xu, B. H.; Zhang, S. J.; Yu, A. B.; Zulli, P. Discrete particle simulation of gas fluidization of particle mixtures. AIChE J. 2004, 50, 1713. (33) Li, Y.; Yang, G. Q.; Zhang, J. P.; Fan, L.-S. Numerical studies of bubble formation dynamics in gas-liquid-olid fluidization at high pressures. Powder Technol. 2001, 116, 246. (34) Xu, B. H.; Yu, A. B.; Chew, S. J.; Zulli, P. Simulation of the gas-solid flow in a bed with lateral gas blasting. Powder Technol. 2000, 109, 14.

(35) Xu, B. H.; Yu, A. B.; Zulli, P. The effect of interparticle forces on powder fluidization behaviour. In Powders and Grains 2001; Kishino, Y. A. A., Ed.; Balkema Publishers: Rotterdam, The Netherlands, 2001; pp 577-580. (36) Xu, B. H.; Feng, Y. Q.; Yu, A. B.; Chew, S. J.; Zulli, P. Particle scale modeling of the gas-solid flow in a fluid bed reactor: Simulation versus experiment. Powder Handl. Process. 2001, 13, 71. (37) Zhang, J.; Li, Y.; Fan, L.-S. Numerical studies of bubble and particle dynamics in a three-phase fluidized bed at elevated pressures. Powder Technol. 2000, 112, 46. (38) Zhang, J.; Li, Y.; Fan, L.-S. Discrete phase simulation of gas-liquid-solid fluidization systems: Single bubble rising behavior. Powder Technol. 2000, 113, 310. (39) Li, J.; Mason, D. J. A computational investigation of transient heat transfer in pneumatic transport of granular particles. Powder Technol. 2000, 112, 273. (40) Ergun, S. Fluid flow through packed columns. Chem. Eng. Prog. 1952, 48, 89. (41) Wen, C. Y.; Yu, Y. H. Mechanics of fluidisation. Chem. Eng. Prog. Symp. Ser. 1966, 62, 100. (42) Di Felice, R. The voidage function for fluid-particle interaction systems. Int. J. Multiphase Flow 1994, 20, 153. (43) Koch, D. L.; Hill, R. J. Inertial effects in suspension and porous-media flows. Annu. Rev. Fluid. Mech. 2001, 33, 619. (44) Xu, B. H.; Yu, A. B. Reply to comments on our paper ‘Numerical simulation of the gas-solid flow in a fluidised bed by combining discrete particle method with computational fluid dynamics’ by Hoomans et al. Chem. Eng. Sci. 1998, 53, 2646. (45) Crowe, C.; Sommerfeld, M.; Tsuji, Y. Multiphase Flows with Droplets and Particles; CRC Press: Boca Raton, FL, 1992. (46) Feng, Y. Q.; Yu, A. B. Comments on ‘Discrete particlecontinuum fluid modelling of gas-solid fluidised beds’ by Kafui et al. (Chem. Eng. Sci. 2002, 57, 2395-2410). Chem. Eng. Sci. 2004, 59, 719. (47) Stokes, G. G. On the effect of the internal friction of fluid on the motion of pendulums. Trans. Cambridge Philos. Soc. 1851, 9, 8. (48) Liang, S.-C.; Hong, T.; Fan, L.-S. Effects of particle arrangements on the drag force of a particle in the intermediate flow regime. Int. J. Multiphase Flow 1996, 22, 285. (49) You, C.; Qi, H.; Xu, X. Drag force in dense gas-particle two-phase flow. ACTA Mech. Sin. 2003, 19, 228. (50) Jasberg, A.; Koponen, A.; Kataja, M.; Timonen, J. Hydrodynamical forces acting on particles in a two-dimensional flow near a solid wall. Comput. Phys. Commun. 2000, 129, 196. (51) Yu, A. B. Powder ProcessingsModels and Simulations (MS 556). In Encyclopedia of Condensed Matter Physics; Elsevier: New York, 2005, in press. (52) Hu, H. H. Direct simulation of flows of solid-liquid mixtures. Int. J. Multiphase Flow 1996, 22, 335. (53) Pan, T.-W.; Joseph, D. D.; Bai, R.; Glowinski, R.; Sarin, V. Fluidization of 1204 spheres: Simulation and experiment. J. Fluid Mech. 2002, 451, 169. (54) Cook, B. K.; Noble, D. R.; Williams, J. R. A direct simulation method for particle-fluid systems. Eng. Comput. 2004, 21, 151. (55) Kafui, K. D.; Thornton, C.; Admas, M. J. Reply to comments by Feng and Yu on ‘Discrete particle-continuum fluid modeling of gas-solid fluidized beds’ by Kafui et al. Chem. Eng. Sci. 2004, 59, 723. (56) Koch, D. L.; Sangani, A. S. Particle pressure and marginal stability limits for a homogeneous monodisperse gas fluidized bed: Kinetic theory and numerical simulations. J. Fluid Mech. 1999, 400, 229. (57) Happel, J. Viscous flow in multi-particle systems: Slow motion of fluids relatice to beds of spherical particles. AIChE J. 1958, 4, 197. (58) Macdonald, I. F.; Ei-Sayed, M. S.; Mow, K.; Dullien, F. A. L. Flow through porous mediasThe Ergun equation revisited. Ind. Eng. Chem. Fundam. 1979, 18, 199. (59) Bouillard, J. X.; Lyczkowski, R. W.; Gidaspow, D. Porosity distributions in a fluidized bed with an immersed obstacle. AIChE J. 1989, 35, 908.

8390


(60) Boemer, A.; Qi, H.; Renz, U. Eulerian simulation of bubble formation at a jet in a two-dimensional fluidized bed. Int. J. Multiphase Flow 1997, 23, 927. (61) van Wachem, B. G. M.; Schouten, J. C.; van den Bleek, C. M.; Krishna, R.; Sinclair, J. L. Comparative analysis of CFD models of dense gas-solid systems. AIChE J. 2001, 47, 1035.

(62) Lacey, P. M. C. Developments in the theory of particle mixing. J. Appl. Chem. 1954, 4, 257.

Received for review July 13, 2004 Revised manuscript received September 16, 2004 Accepted September 30, 2004 IE049387V

Assessment of Model Formulations in the Discrete Particle Simulation

Recommend Documents