Chemical Engineering Science 91 (2013) 1–4
Contents lists available at SciVerse ScienceDirect
Chemical Engineering Science journal homepage: www.elsevier.com/locate/ces
Note
Fully resolved simulation of a gas-fluidized bed: A critical test of DEM models S.H.L. Kriebitzsch n, M.A. van der Hoef, J.A.M. Kuipers Department of Chemical Engineering & Chemistry, Eindhoven University of Technology, P.O. Box 513, 5600MB Eindhoven, The Netherlands
H I G H L I G H T S c c c
Direct Numerical Simulation of gas-fluidized beds were done. The individual fluid–particle drag is compared with the force as evaluated in unresolved Discrete Element Models. The average DEM drag is found to be about 33% smaller than the ‘‘true’’ DNS drag.
a r t i c l e i n f o
abstract
Article history: Received 14 May 2012 Received in revised form 6 December 2012 Accepted 20 December 2012 Available online 31 December 2012
We performed Direct Numerical Simulation (DNS) of a gas-fluidized bed using the Immersed Boundary method combined with traditional Computational Fluid Dynamics. We have compared the individual fluid–particle force with the force as evaluated in an unresolved Discrete Element Model (DEM) simulation using closures for the gas–solid force for one fluidisation condition. We find that the average DEM gas–solid force is about 33% smaller than the ‘‘true’’ value which follows from the DNS model. For more realistic DEM simulations, a modification of the gas–solid force calculation is required, for instance by adding fluctuations or including the effect of the granular temperature. & 2012 Elsevier Ltd. All rights reserved.
Keywords: Fluidization Multiphase flow Particle Fluid mechanics Direct numerical simulation Gas–solid force
1. Introduction In the past 10 years, discrete elements models (DEMs) have become increasingly popular for modeling gas fluidized beds, and have now established themselves as a standard tool in the study of such systems (Deen et al., 2007) on a laboratory scale. A typical feature of discrete element models is that the flow is unresolved, that is, the gas phase is described on a scale larger than the size of the particles. As a consequence, the gas–solid interaction cannot be described at the scale of a single particle, but rather at the scale of the smallest computational unit of the gas-phase, i.e. a single computational fluid dynamics (CFD) cell. The gas–solid interaction force on a particle i with velocity vi is then determined by the cellaveraged porosity e and flow velocity u, typically via a relation like FDEM g,i ¼ 3pmdpðuvi Þ Fðe,Rei Þ,
Re ¼
rde9uvi 9 , m
n
Fg ¼ Fd V P =p,
ð2Þ
where VP is the volume of the particle. We consider two correlations for Fðe,ReÞ: the combined Ergun and Wen and Yu correlation (see Gidaspow, 1994):
ð1Þ
Corresponding author. Tel.: þ31 40 247 3666; fax: þ31 40 247 5833. E-mail addresses:
[email protected] (S.H.L. Kriebitzsch),
[email protected] (M.A. van der Hoef),
[email protected] (J.A.M. Kuipers). 0009-2509/$ - see front matter & 2012 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.ces.2012.12.038
where F is a correlation that has to specified, d is the particle diameter (for simplicity we consider monodisperse systems), and m and r are the viscosity and density of the gas phase, respectively. Note that it is common to evaluate the porosity and flow velocity at the location of particle, by an interpolation of the values of e and u of all the cells within a certain range of the particle, so that the values depend slightly on the location of the particle (Deen et al., 2007). We stress that Eq. (1) is the total gas–solid interaction force, which can be split into a drag force and a pressure gradient force:
Fðe,ReÞ ¼
8 3:65 > > 1 þ0:15Re0:687
for e 4 0:8,
150 ð1eÞ 1:75 Re > > þ : 18 e2 18 e2
for e o 0:8
ð3Þ
2
S.H.L. Kriebitzsch et al. / Chemical Engineering Science 91 (2013) 1–4
and the correlation by Beetstra et al. (2007), obtained from lattice Boltzmann simulations (with f ¼ 1eÞ: " # 180 f 0:413 e1 þ 3ef þ 8:4Re0:343 Re 1=2 2 Fðe,ReÞ ¼ þ e ð1 þ 1:5f Þ þ : 18 e2 24 e2 1þ 103f Reð1 þ 4fÞ=2 ð4Þ In this paper we want to put the current class of DEMs to the test, by performing a fully resolved simulation of a fluidized bed – which treats the gas–particle interaction ab initio. Since in such a simulation the ‘‘true’’ force on a particle is known at any moment in time, we can compare this with the gas–solid force that would have been evaluated for that particle in the framework of a DEM simulation using relations of the type (1). To our knowledge, this is the first fully resolved simulation of an actual gas-fluidized bed. 2. Simulation details 2.1. DNS method For the solid phase, a simple hard sphere code is used, where the particle collisions are binary and instantaneous; the change in velocity and angular momentum due to the gas–solid force is calculated every CFD timestep. Note that this is different from DEM simulations, where only the change in linear momentum due to the gas–solid interaction is modelled. The gas flow field is obtained from the full Navier–Stokes equation using state-of-the art Computational Fluid Dynamics (CFD) methods, where a firstorder-accurate method is used to discretize the momentum equations in time, in which the pressure, viscous and convective terms are treated implicit, semi-implicit and explicit, respectively. The noslip boundary conditions at the surface of the particle are enforced by use of the Immersed Boundary (IB) method (Peskin, 2002); This is IB done by adding a force term f to the Navier–Stokes equations: @ IB u þ = uu ¼ =P þ mDuþ f , rG ð5Þ @t with fluid density rG , fluid velocity u, pressure P and dynamic viscosity m. The implementation that we adopt is along the lines of Uhlmann (2005). We have validated the code for a number of standard test cases, such as the hydrodynamics interaction force between two approaching spheres, and the hydrodynamic force for flow past regular and random arrays. Here only a brief outline of the method is given, details of the implementation and verification are published in Kriebitzsch (2011). The basic idea is that the surface of the sphere is covered with markerpoints (see Fig. 1), each of which exerts a force on the fluid such that the gas-phase velocity is equal to the surface velocity at the location of the marker point, thereby modeling ‘‘no-slip’’ boundary IB conditions. The ‘‘essential’’ steps to compute the IB forcing term f are as follows:
interpolation of fluid velocity to markerpoints m: X Um ¼
Dðxi,j,k X m Þ ui,j,k
ð6Þ
i,j,k
compute the force at the markerpoint as the difference of the current interpolated fluid velocity and the velocity of the particle at the location of the markerpoint: F IB m¼
r Dt
ðV m U m Þ
ð7Þ
mapping of force F IB back to the Eulerian grid: IB
f i,j,k ¼
X DV m Dðxi,j,k X m Þ F IB m b m
d
Fig. 1. The 2D sketch of an immersed object. The object is represented by a IB number of markerpoints, which are used to compute the force density f .
Dðxi,j,k X m Þ denotes a so-called discrete or regularized delta function with a finite support, which is computed as the product of one-dimensional regularized delta function. As indicated by the shaded squares in Fig. 1, we use a width of three Eulerian grid cells as support, hence in a 3 3 3 cube in 3D. The sum of the forces F IB from the markerpoint of one sphere i is equal to the total force that the sphere applies to the fluid, and hence the negative of this is FDNS g,i , which is used in the update of the position of the sphere. 2.2. Simulation settings The physical properties of air are used for the density and dynamic viscosity of the gas phase: r ¼ 1:2 kg=m3 , m ¼ 1:8 105 Pa.s. The particle diameter d and particle density rp are set to d¼1 mm and rp ¼ 1000 kg=m3 , which in the Geldart classification this corresponds to B-type particles. The timestep to update the gas flow was set to dt ¼ 1 105 s, while the gridsize dl was set to d=10. Choosing a coarser grid (i.e. d/6.5) did not influence the results. We performed simulation of a gas-fluidized bed of 2000 particles, with a width equal to W¼14.8, d¼14.8 mm, a depth D¼5.6, d¼5.6 mm, and a height H¼42 mm; the height of the bed in the packed state was 21.2 mm. The bed was fluidized with a superficial velocity U o ¼ 0:5 m=s, which is roughly two times the minimum fluidization velocity. The typical particle Reynolds number at this velocity is Re ¼ rdU o =m ¼ 35:3. Constant pressure outflow boundary condition was used at the top, and no-slip boundary conditions were used on all side walls. In the simulations, relatively homogeneous fluidisation was observed with a wave-like expansion behaviour of the particle bed that means the bedheight fluctuates periodically (see also Verloop and Heertjes, 1974). The average bed height in the fluidized state was 29.0 mm, and the corresponding porosity was 0.563. For comparison, we have also run a DEM simulation of the same system; the timestep was the same as in the DNS and the grid size was about 3d. We found that the bedheight was significantly lower in that case (26.2 mm), and also the characteristics of the bed oscillations was different. Snapshots from both simulations are shown in Fig. 2. 3. Comparisons with DEM type gas–solid force models
ð8Þ
One of the assets of DNSs is that the actual (or ‘‘true’’) gas– solid force on each individual particle is known at any point in
S.H.L. Kriebitzsch et al. / Chemical Engineering Science 91 (2013) 1–4
3
Fig. 2. Snapshots from the DNS (middle) and DEM (right) simulation of a fluidized bed with NP ¼ 2000 particles. The fluid velocity vectors (see enlarged section of the DNS on the left) very clearly bring out the differences in the detail of modelling in the DNS and DEM simulations.
time, which allows us to address the key question of this work: how well do DEM-type gas–solid force models predict the actual force acting on a particle? Or in other words: what is the deviation of the actual force FDNS acting on particle i as computed g,i in a fully resolved simulation, from the force FDEM one would g,i compute from a relation like (1)? Very recently this was investigated for an ideal, well-defined system: static particles, periodic boundary conditions, homogeneous arrays (Kriebitzsch et al., 2013); in this work the same is done for an actual fluidized bed system. Note that for the comparison that we make here, we do not compare DNS results with the results of a DEM simulation; rather, at any point in time of the DNS, a gas–solid force is computed for each particle as if it were a DEM simulation, that is, by using relation (1). For this, the detailed scale information on the flow field of the DNS has been transformed to a coarse scale flow field and porosity field as available in a DEM-type simulation. We have set the grid size of this coarser ‘‘DEM grid’’ to 3d. From the DNS we have evaluated histograms for DF i ¼ DNS DNS 9FDEM and FDEM g,i 99Fg,i 9 from N ¼120 000 ‘‘samples’’ for Fg,i g,i , once the fluidisation had reached steady-state. The samples are the data for Fg,i for each of the 2000 particles at 60 different time steps, each separated by an interval of 1000 time steps. In Fig. 3 (left graph) we show the histogram for DF i normalized by the magnitude of the buoyancy corrected gravity force 3 9F g F b 9 ¼ 16pd ðrp rÞg, where the DEM gas–solid force is evaluated both from the correlation (4) by Beetstra et al. (2007) (indicated by BHK), and from correlation (3) by Ergun and Wen and Yu. For most particles the computed DEM-type fluid–particle force is substantially lower than the ‘‘true’’ force obtained from the DNS. Since the system is in the fluidized state, the average
force /9FDNS 9S is about equal to the buoyant ‘‘weight’’ F g F b of the particles, however, the average value of the computed DEMtype fluid–particle force is about 33% lower. We find also that the DEM gas–solid force calculated by (4) and (3) gives almost identical results. This is not surprising, since at this particular average bed porosity and Reynolds number (e ¼ 0:563, Re ¼35.3) the correlations (3) and (4) given almost identical values for F: 39.64 and 40.34. For other Reynolds numbers we expect a more pronounced difference: at the same porosity (e ¼ 0:563), the Ergun/Wen and Yu force is 20% smaller for Re ¼1, and 20% larger for Re ¼125, compared to the Beetstra–vdHoef–Kuipers force (4). In the comparison discussed above we have used the individual particle velocity vi in the evaluation of FDEM (see Eq. (1)). One g,i could argue that using the cell-average solid velocity v instead of vi is more consistent with how the fluid–particle force correlations have been derived, since correlations like (3) and (4) represent the hydrodynamic force for a steady state; the instantaneous force on particle i will fluctuate with fluctuating vi , but these fluctuations will not necessarily scale in the same way as the average hydrodynamic force scales with the average slip velocity. In Fig. 3 the triangles represent the difference of the DNS fluid–particle force with the DEM fluid–particle force where P cell we used v ¼ ð1=N cell Þ N v for each particle in the DEM cell, i A cell i instead of vi . Interestingly, we see that the effect on the DEM fluid–particle force is negligible.
4. Discussion Our main finding is that on average the DEM force is significantly lower than the (true) DNS force, consistent with the
4
S.H.L. Kriebitzsch et al. / Chemical Engineering Science 91 (2013) 1–4
1 BHK, vi BHK, v
0.10 0.1
0.05
0.00 -1.5
-1
-0.5 0 ΔFi / (Fg-Fb)
0.5
1
0
0.2
0.4
δF
0.6
fraction of particles with δFi > δF
Ni / N
min
Ergun, vi
0.15
0.01 0.8
min
DNS Fig. 3. Left: Histogram of the difference DF i ¼ 9FDEM g,i 99Fg,i 9, normalized by the buoyant weight of one particle, computed from a fully resolved simulation of a fluidized bed with 2000 particles. N i =N is the fraction of particles that have a corresponding specific value of DF i . N¼ 20 000 particles were used to compute the histogram by taking DEM,scaled min 99FDNS , sampling from 60 different points in time. Right: Fraction of particles for which the relative deviation dF i ¼ ð9Fg,i g,i 9Þ=ðF g F b Þ is – on average – larger than dF DEM,scaled where Fg,i is the DEM force scaled such that its average is equal to the average DNS force.
reduced bedheight observed in the direct comparison between the DEM simulation and DNS. Considering that typical DEM-type closures for the hydrodynamic force are obtained from simulations or experiments with beds of stationary particles, this suggests that the mean gas–solid force in systems with a granular temperature, such as fluidized beds, is higher than for a particle bed with the same porosity and stationary particles. An increase of the dimensionless gas–solid force with the granular temperature have also been reported by Whylie et al. (2003), although for a different kind of systems with a porosity of e ¼ 0:8. Our simulations also showed that the effect of the granular temperature is not captured by using the individual particle velocity vi in the gas–solid force correlation, at least not for this particular system. It should be stressed that these conclusions are based on the results for only one fluidization condition, since the DNSs are currently very time consuming. Clearly far more different systems should be studied in order to make more definite and quantitative statements, and validated with dedicated one-to-one experiments. However, from the present study it is already clear that accurate correlations for static, random systems are not a priori applicable to fluidized systems, and that one can make no simple statement on whether the Ergun and Wen and Yu or the Beetstra– vdHoef–Kuipers closure for the gas–solid force is better suited for such systems. For instance, for Re 450, the Ergun gas–solid force is higher than the gas–solid force computed with Beetstra– vdHoef–Kuipers (Eq. (4)), and hence would be closer to the DNS result. For Re o 10, the reverse is true, and the Beetstra–vdHoef– Kuipers force will be more accurate. From this it is clear that the appropriateness of these static correlations for fluidized systems is a bit of a hit and miss affair. An engineering approach would be to perform DNSs for a wide range of fluidized beds, measure the average gas–solid force, and determine the scaling required for the closures given by Ergun or Beetstra–vdHoef–Kuipers accordingly; the drawback is that it would provide a correlation that is not based on physical principles, and its use limited to fluidized systems. A more fundamental approach would be to assume that
one of the key parameter which modifies the gas–solid force for fluidized systems is the local granular temperature y; so the approach would be to study the effect of y in detail by performing DNSs for small periodic systems with a well-defined granular temperature, in the spirit of the closures derived for static systems. We expect that a combination of both approaches would provide a major step forward in the discrete element modeling of gas–fluidized systems. Apart from the fact that the mean force is quite different for the DEM and DNS model, there is also the issue of the rather large scatter of the DNS force; that is, even when the average DEM gas– solid force would be equal to the DNS result (via scaling or other methods), still for one out of three particles the individual DEM gas–solid force deviates more than 25% with the true individual DNS gas–solid force, and for one out of 10 particles the deviation is more than 40% (see Fig. 3, right graph). This spreading is due to temporal and spatial fluctuations in the particle positions and velocities which are, by construction, not captured in DEM.
References Beetstra, R., van der Hoef, M.A., Kuipers, J.A.M., 2007. Drag force of intermediate reynolds number flow past mono- and bidisperse arrays of spheres. AIChE J. l 53, 489–501. Deen, N.G., van Sint Annaland, M., van der Hoef, M.A., Kuipers, J.A.M., 2007. Review of discrete particle modeling of fluidized beds. Chem. Eng. Sci. 62, 28–44. Gidaspow, D., 1994. Multiphase Flow and Fluidization. Academic Press, San Diego. Kriebitzsch, S.H.L., 2011. Direct Numerical Simulation of Dense Gas-Solid Flows. Ph.D. Thesis. Eindhoven University of Technology. Kriebitzsch, S.H.L., van der Hoef, M.A., Kuipers, J.A.M., 2013. Drag force in discrete particle models - continuum scale of single particle scale? AIChE J. 59, 316–324. Peskin, C.S., 2002. The immersed boundary method. Acta Numer. 11, 479–517. Uhlmann, M., 2005. An immersed boundary method with direct forcing for the simulation of particulate flows. J. Comput. Phys. 209, 448–476. Verloop, J., Heertjes, P.M., 1974. Periodic pressure fluctuations in fluidized beds. Chem. Eng. Sci. 29, 1035–1042. Whylie, J.J., Koch, D.L., Ladd, A.J.C., 2003. Rheology of suspensions with high particle inertia and moderate fluid inertia. J. Fluid Mech. 480, 95–118.