Computers & Fluids 68 (2012) 159–167
Contents lists available at SciVerse ScienceDirect
Computers & Fluids j o u r n a l h o m e p a g e : w w w . e l s e v i e r . c o m / l o c a t e / c o m p fl u i d
Numerical simulation of heat and mass transfer in pneumatic conveying dryer Samy M. El-Behery a,⇑, W.A. El-Askary a, Mofreh H. Hamed b,c, K.A. Ibrahim a a
Faculty of Engineering, Menoufiya University, Shebin El-kom, Egypt Faculty of Engineering, Kafrelsheikh University, Kafrelsheikh, Egypt c Faculty of Engineering, Islamic University in Madina, Madina, Saudi Arabia b
a r t i c l e
i n f o
Article history: Received 29 December 2011 Received in revised form 30 June 2012 Accepted 11 August 2012 Available online 20 August 2012 Keywords: Heat and mass transfer Pneumatic dryer Two-phase Numerical simulation
a b s t r a c t Two-dimensional Eulerian–Lagrangian model is presented for heat and mass transfer in pneumatic conveying dryer. The model takes into account the particle–particle and particle–wall collisions, lift forces, particle rotation, turbulence modulation and turbulence dispersion (i.e., four-way coupling). The drying simulation is based on a two-stage drying model. Different correlations for heat transfer coefficient are tested and assessed in terms of their accuracy. The model is validated against the available experimental data and good agreement is obtained. The model predictions are compared to other models from literature and it produces better results than existing models. It is also found that the turbulence dispersion has greater effect on the model predictions than particle–particle collision. However, neglecting either particle–particle collision or turbulence dispersion results in a lower heat transfer and drying rates. Ó 2012 Elsevier Ltd. All rights reserved.
1. Introduction Drying is an essential process in many industrial applications such as food, agricultural, ceramic, polymers and plastic, pulp and paper, pharmaceutical and wood processing industries. Drying equipments are classified according to the heat transfer mechanism to convective and conductive and according to the handling characteristics to batch and continuous operation [1]. Pneumatic conveying dryers which can be classified as convective and continuous drying equipment is one of the most widely used equipment. Pneumatic dryers are characterized by simultaneous momentum, heat and mass transfer processes between the dispersed material and the drying agent. The large surface area for heat and mass transfer results in higher drying rate and higher drying capacity. In these types of dryers the contact time between the drying medium and particulate material is relatively short (usually few seconds only). Therefore, these dryers are suitable for heat-sensitive materials and also for removing external moisture. This allows higher inlet temperatures to be used than in many other dryers without unduly heating the product [1]. Pneumatic dryers are simple in construction and have low capital cost. Vertical type of construction, which facilitates installation in exiting buildings, is an advantage of pneumatic dryer systems [2]. Among other dryers, the pneumatic dryer shows the highest removal rate of the liquid from the solid particle [3]. Drying calculations are based on the knowledge of air and material properties. The successful ⇑ Corresponding author. E-mail address:
[email protected] (S.M. El-Behery). 0045-7930/$ - see front matter Ó 2012 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.compfluid.2012.08.006
simulation of a pneumatic dryer requires knowledge of the solids residence time (and hence the particle motion within the dryer), the heat transfer coefficient and mass transfer process inside and outside the solid. By employing a volumetric heat transfer concept, as used for rotary dryers, simple estimation procedures have been suggested by Williams-Gardner [4] and Moyers and Baldwin [5]. These procedures assumed that the particles are traveling at a steady velocity close to the gas velocity. Baeyens et al. [1] pointed out that these methods can over-predict the required dryer length by 200–400%. To model the acceleration zone near the feed point a stepwise procedure was suggested by many researchers such, Thorpe et al. [6], Kemp et al. [7], Radford [8] and Kemp and Oakley [9]. Although these procedures are considerable improvement on the steadystate, Kemp et al. [10] reported that they still give errors of 50–100% in the tube length prediction. Further development in the dryer simulations is the application of one-dimensional model based on the two-fluid theory. Levy and Borde [11], Pelegrina and Crapiste [12] and El-Behery et al. [13] simulated the pneumatic dryer using the one-dimensional steady-state model. However, this model does not provide any information about the cross section distribution of flow properties, which may be an important aspect to be considered in drying processes [14]. To predict the cross-sectional distribution of flow properties Skuratovsky et al. [15,16] developed a two-dimensional steady-state model based on the two-fluid theory. The predictions of the two-dimensional model did not present any significant difference as compared to those provided by Levy and Borde [11]. Recently, the computational fluid dynamics (CFD) code, Fluent, was used to perform 3-D simulation of pneumatic
160
S.M. El-Behery et al. / Computers & Fluids 68 (2012) 159–167
Nomenclature Symbols Cp Cl D Dv dp FD FLR FSL Hfg h hh hm Ip M o m m Nu Pr Prt Rt r T To t
specific heat (J/kg K) constant in turbulence model (–) pipe diameter (m) diffusion coefficient (m2/s) particle diameter (lm) drag force (N) Magnus lift force (N) shear lift force (N) latent heat of vaporization (J/kg) enthalpy (J/kg) heat transfer coefficient (W/m2 K) mass transfer coefficient (m/s) particle moment of inertia (N m) molecular weight (kg/kmole) mass flow rate (kg/s) mole fraction (–) Nusselt number (–) Prandtl number (–) turbulent Prandtl number (–) turbulent Reynolds number (–) radial coordinate (m) temperature (K) torque (N m) time (s)
drying process [17–19]. Eulerian–Lagrangian approach was applied by Mezhericher et al. [17] while Eulerian–Eulerian approach was used by Jamaleddine and Ray [18,19]. The model of Mezhericher et al. [17] can not predict correctly either the gas temperature or particle water content. On the other hand, the predictions of Jamaleddine and Ray [18,19] agree well with the experimental data. This can be attributed to the limitations of the available Eulerian– Lagrangian model in Fluent. Experimental investigations of the pneumatic dryer were given also by NamKung and Cho [20] and Kaensup et al. [21,22]. A critical survey of the earlier work reveals a good amount of theoretical work on pneumatic conveying dryer. However, the work devoted to computational fluid dynamics including turbulence modeling is very scarce. In addition, the effects of many important phenomena such as particle–particle collision, particle–turbulence interaction, lift forces, particle size distribution and turbulence modulation are limited. Chagras et al. [23], Boulet and Moissette [24], Mansoori et al. [25] and El-Behery et al. [26] showed that these parameters have a great effect on the gas–solid heat transfer. Furthermore, experimental work on pneumatic conveying dryer is limited. Therefore, the main objective of this paper is to provide a validated computational model, including those parameters that previously neglected, that can be used for pneumatic conveying design and operation predictions.
X x Y
solid water content (kgwater/kgsolid) axial coordinate (m) water vapor mass fraction (kgwater/kgair)
Greek symbols a gas void fraction (–) b solids void fraction (–) d porosity (–) k thermal conductivity (W/m K) l viscosity (kg/m s) lt turbulent viscosity (kg/m s) q density (kg/m3) R universal gas constant (J/kmole K) sw shear stress at wall (N/m2) x angular velocity (1/s) Subscripts Ave average cr critical da dry air g gas H2O water vapor p particle v vapor
2.1. Model assumptions 1. The flow model is for a 2-D/axisymmetric duct. 2. The particles are spherical in the particle tracking procedure. 3. The particles are rigid during particle–particle and particle–wall collisions. 4. The gas phase is assumed as an ideal gas-mixture of air and water vapor. 5. Heat and mass transfer occurs between gas phase and individual particles (i.e., heat and mass transfer between particles themselves is ignored). The particle-to-particle heat transfer is neglected based on the finding of Mansoori et al. [27]. They reported that the particle-toparticle heat transfer can result from three main mechanisms, namely: heat transfer by radiation, heat transfer through the contact point, and heat exchange through the gas lens at the interface between colliding particles. They reported also that the first mechanism is only significant at temperatures higher than 600 °C. They also found in their numerical efforts, using hot and cold particles streams, that the inter-particles heat transfer does not significantly affect the mean gas and particles temperature, and becomes important only when the rate of change in temperature of hot or cold particles is important. 2.2. Gas flow modeling
2. Mathematical model Four-way coupling Eulerian–Lagrangian approach is used in the present study to predict the heat and mass transfer in gas–solid flow through pneumatic conveying dryer. The gas phase is simulated using Reynolds Averaged Navier–Stokes Equations. The turbulent viscosity is modeled by the low-Re k–e model and the particle tracking procedure is used for the solid phase. To provide a reasonable solution for engineering objectives some simplifying assumptions are taken as follow.
The general form of the elliptic differential equations governing axisymmetric, turbulent, steady, compressible and non-isothermal two-phase flow is:
@ @/ 1 @ @/ þ ¼ S/ S/p aqu/ aC/ aqrv / aC/ r @x @x r @r @r
ð1Þ
where S/ and S/p are source terms of gas and dispersed phases respectively, while exchange coefficient, C/ is summarized in Table 1 for the dependent variable /.
161
S.M. El-Behery et al. / Computers & Fluids 68 (2012) 159–167 Table 1 Governing equations of the gas phase. Conservation of
/
C/
S/
Mass Axial momentum
1 Ug
0.0
0.0
Radial momentum
Vg
leff
Energy equation
Tg
lg
Turbulent kinetic energy
k
Dissipation rate
e
leff rk leff re
Water vapor mass fraction
Y H2 O
Pr
2
k
ð2Þ
e
In the low-Re k–e model damping functions and the extra source terms D and E (see Table 1) are used instead of the wall function. These source terms and damping functions are calculated as given by Launder and Sharma [28].
f 2 ¼ 1 0:3 exp R2t
f1 ¼ 1:0;
D ¼ 2l
pffiffiffi!2 @ k ; @y
E¼2
llt @U g q @y
ð3Þ 2 ð4Þ
To account for the increase of turbulent Prandtl number near the wall, the expression given by Kays [29] is used as follows:
1 l ¼ 0:5882 þ 0:228 t 0:0441 Prt lg lg 1 exp 5:165
lt lg
!2
ð5Þ
lt
2.3. Particle phase modeling The solid phase is simulated using the Lagrangian approach, a few thousands of computational particles ‘parcels’ are traced through the flowfield in each coupling iteration. After each given time step the new position of the parcels and the new transitional and angular velocities are calculated from the equations of motion as in [26] through:
! dxp ! ¼ up dt
Ip
aðG qe DÞ a ke ðC e1 f1 G C e2 f2 qeÞ þ aE 0.0
nus lift due to particle rotation and gravity, respectively. Detailed description of these forces can be found in [26]. The particle temperature is traced along the particle trajectory by the following ordinary differential equation:
mp C pp
o dT p 2 ¼ hh pdp ðT g T p Þ md Hfg dt
ð10Þ
In order to solve the forgoing equations (i.e. Eqs. (6)–(10)); the instantaneous gas velocity components at the particle’s location are needed. The mean velocities are interpolated from neighboring grid points, while the fluctuating components are generated using the Stochastic Separated Flow model (SSF) given by Shuen et al. [30]. 2.4. Mass transfer model The mass transfer in the present model is based on the twostage drying process. In the first stage, the solid surface can be considered to be fully wetted and the resistance to the mass transfer is located in the gas side. The evaporation rate from individual particle can be expressed as given by Levy and Borde [11] as: o
The turbulence model constants are, Ce1 = 1.44, Ce2 = 1.92, Cl = 0.09, rk = 1.0 and re = 1.3.
mp
l
þ Prtt
lt qDv þ Sc t
The effective and eddy viscosities are calculated as:
leff ¼ l þ lt ; lt ¼ C l fl q
h i ap g @ 2 @ ~ @@x þ @x aleff @U ar leff @V@xg þ 1r @r @x 3 r V g h i ap g @ @ þ 1r @r @@r þ @x aleff @U aleff ðr @V@rg 23 r ~ VgÞ @r V Vg 2aleff r2g þ 23 alr e r ~ h i 2 1 @ @ 2C pg @x aqg U g ðU g þ V 2g Þ þ 1r @r ðr aqg V g ðU 2g þ V 2g ÞÞ
leff
2
md ¼ hm pdp
MH2 O pv o MH2 O pv g RT p RT g
ð11Þ
where pvo and pvg are the partial pressures of water vapor at the particle surface and the gas phase. Mezhericher et al. [17] introduced the following differential equation to calculate the time change of particle radius during the constant rate period. o dRp 1 md ¼ dt qw 4pR2p
ð12Þ
ð6Þ
! ! ! ! dup ! ¼ F D þ F LS þ F LR þ F g dt
ð7Þ
! dxp ! ¼ T dt
ð8Þ
! ! 3 ! To ¼ pldp ½0:5r U g x p
ð9Þ
! ! ! where x p is the particle position vector, U g ; u p are the gas and par! ticle velocity vectors, x p is the particle angular velocity vector, To is 2
the torque acting on the particle, Ip ¼ 0:1mp dp is the particle mo! ! ! ! ment of inertia and mp is the particle mass. F D ; F LS ; F LR and F g are the components of the forces arising from drag, shear lift, Mag-
(a) Constant rate period
(b) Falling rate period
Fig. 1. Two-stage drying model of wet particle.
162
S.M. El-Behery et al. / Computers & Fluids 68 (2012) 159–167
In the present study an analytical expression is developed to calculate the instantaneous particle diameter as a function of dry particle diameter do, and the instantaneous solid water content, X, as follows:
dp ¼ do
Frantz correlation [9].
1=3
qs ðX X cr Þ þ 1 qw
ð13Þ
The second drying stage period starts when the particulate surface becomes no longer wetted and evaporation must occur from within the pores and dry porous crust starts formation around the wet core, as shown in Fig. 1. This is assumed to occur at solid water content, X, less than the critical solid water content, Xcr. The particle diameter, do will no longer changed during this stage while the radius of the wet core recedes to the particle center. Heat is conducted from the particle surface through the dry crust to the interface and water evaporates and diffuses back to the particle surface. The rate of vapor diffusion from the wet core to the particle surface decreases as the dry crust thickens increases. This process is called falling rate period. By introducing the resistance in the dry porous crust given by Maxwell–Stefan type diffusion, the mass transfer rate can be calculated as given by Skuratovsky et al. [15,16] as: o
md ¼
1
1 1 4pd Ro Ri
1
M H2 O D v P P P v i ln P Pv o RT av e
ð14Þ
The mass transfer rate from the wet core, Eq. (14), must be equal to that from particle surface to the gas, Eq. (11). Eliminating the partial pressure of water vapor at the dry crust, pvo from Eq. (14) using Eq. (11) yields an implicit equation for the evaporation rate: o
md ¼
1 1 1 4pd Ro Ri
ln P
1
M H2 O D v P RT av e P Pv i o
RT p MH
h 2O m
d2po
p
T
md þ pv g T pg
ð15Þ
The heat to the wet core is transferred first by convection from the ambient gas and then by conduction through the dry crust. Assuming quasi-steady state conditions, the temperature of the wet core can be calculated by:
ðT g T p Þ
!1
1 2
hh pdo
¼ Tg Ti
1
do di þ 2 hpdo 2pdo di ks
!1 ð16Þ
where Tp, Ti are temperatures of the outer dry crust and the inner wet core, respectively, and Tave = (Tp + Ti)/2. In general, during the second drying period, the outer shape of the particle might be changed due to shrinkage of both outer and core diameters. However, to simplify the model, it is assumed that the particle outer diameter remains constant during the second drying period. Thus, only the change of the wet core diameter, Di is considered.
dDi 2 o md ¼ dt dpD2i
ð17Þ
2.5. Heat transfer coefficient The convective heat transfer coefficient, hh, was calculated from Nusselt number, Nu, which is expressed as a function of Reynolds number, Rep and Prandtl number, Pr, which are defined as:
Nu ¼
hdp ; kg
Rep ¼
! qg dp j ! ug u p j
lg
Various empirical correlations can be used to calculate the heat transfer coefficient. The following published correlations have been tested in the present study:
;
Pr ¼
l g C pg kg
0:667 Nu ¼ 0:015Re1:6 p Pr
De Brandt correlation [1,11]. 0:667 Nu ¼ 0:16Re1:3 p Pr
ð20Þ
Debrand correlation [31].
Nu ¼ 0:035Re1:15 Pr0:333 p
ð21Þ
Baeyens et al. correlation [1]. The correlation was developed for a large scale pneumatic dryer.
Nu ¼ 0:15Rep
ð22Þ
Modified Ranz–Marshall correlation [11]. The correlation was developed for single droplet/wet particle and it takes into account the resistance of the liquid vapor around the particle to the heat transfer by Spalding number, B.
Nu ¼
B¼
0:333 2 þ 0:6Re0:5 p Pr
ð1 þ BÞ0:7
C pH2 O ðT g T p Þ Hfg
ð23Þ
ð24Þ
Modified Weber correlation [7] An additional term proportional to Re0:8 was added to Ranz– p Marshall correlation to account for turbulent flow.
0:8 Nu ¼ 2 þ 0:5Re0:5 Pr0:333 p þ 0:06Rep
ð25Þ
2.6. Mass transfer coefficient In analogy to the heat transfer coefficient, hh, the mass transfer coefficient, hm, is calculated from Sherwood number, Sh, which is equivalent to Nusselt number, Nu. It is often expressed as a function of the particle Reynolds number, Rep, and Schmidt number, Sc, which is equivalent to Prandtl number, Pr, and they are defined as:
Sh ¼
hm dp ; Dv
Sc ¼
lg qg Dv
ð26Þ
Eqs. (19)–(25) have been tested and used to calculate the mass transfer coefficient in the present study. 2.7. Coupling between the two phases The particles occupy the computational cell and reduce the gas volume fraction. They also exert interaction forces on the surrounding gas phase. Thus, the two phases are coupled through the gas volume fraction and through the total source/sink term, S/P that accounts for the momentum, heat and mass exchange between continuous and dispersed phases. The void fractions of dispersed phase, b and gas phase, a are calculated respectively using trajectory method as depicted in [32] by:
b¼
X N k Dt k V k traj
ð18Þ
ð19Þ
V Cell
;
a¼1b
ð27Þ
where Nk is the number of actual particles in the computational particle ‘parcel’ (k), Vk is the volume of the particle, Vcell is the volume of
163
S.M. El-Behery et al. / Computers & Fluids 68 (2012) 159–167
P
computational cell and traj means summing over all trajectory passing through the computational cell. The source term of dispersed phase in the gas momentum equation is calculated as in [33] by:
SUp i
1 X ¼ mk N k V Cell traj " ! # Nt X q nþ1 n ððupi Þk ðupi Þk Þ g i 1 Dt L
qp
n¼1
ð28Þ
where DtL is the Lagrangian time step used in the solution of Eqs. (3)–(8), and summing over n indicates averaging along particle trajectory (time average). The energy source term, STS , which represents the convective heat exchange between the dispersed phase and the continuous phase as well as the energy transferred to the gas phase by water vapor, is given by:
STp ¼
Nt h i X o 1 X 2 Nk ðhh pdp ðT g T p Þ md hH2 O ÞDtL C pg V Cell traj n¼1
ð29Þ
Smass S
2.9. Supplementary equations In order to solve the above set of equations several supplementary equations, definitions and empirical correlations are required. These will be presented subsequently. It should be noted that both the gas and solid phases are mixtures and hence their thermodynamic properties are calculated using the mixture theory. – Mole fraction of water vapor in the gas stream.
mH2 O ¼
P Rg T g
ð35Þ
R ; Mg
M g ¼ mH2 O M H2 O þ ð1 mH2 O ÞMda
ð36Þ
– Heat capacity of the gas stream.
The effect of particles on the turbulence of the continuous phase is known as turbulence modulation. El-Behery et al. [26] compared the performance of four commonly used turbulence modulation models and they found that the model of Lain and Sommerfeld [33] gives the best agreement with experimental data. This model has the following form and is adopted in the present study:
Cpg ¼ Y H2 O C pH2 O þ ð1 Y H2 O ÞCpda
Skp ¼ upi SUp i U gi SUp i
kg ¼ Y H2 O kH2 O þ ð1 Y H2 O Þkda
ð31Þ
ð34Þ
where
Rg ¼ ð30Þ
Y H2 O Y H2 O þ ð1 Y H2 O ÞM H2 O =M da
– Density of gas stream. The mole fraction of water vapor is used together with the ideal gas equation to calculate the density of the gas phase as follows:
qg ¼
where hH2 O is the enthalpy of water vapor at particle temperature. The mass source term (i.e., the source term in continuity and water vapor mass fraction equations) is calculated as: Nt X o 1 X ¼ Nk ½md Dt L V C traj n¼1
assumed to be 0.25 for particle–particle collision and 0.2 for particle wall-collision.
ð37Þ
– Viscosity of gas stream.
lg ¼ mH2 O lH2 O þ ð1 mH2 O Þlda
ð38Þ
– Thermal conductivity of gas stream.
ð39Þ
– Heat capacity of the dispersed phase. e
Sp ¼ C e3
e k
Skp
ð32Þ
where Ce3 is a model constant, its value varies in the literature from 1.1 to 2. Zhang and Reese [34] accounted for the reduction in the gas turbulence length scale due to the presence of a second phase in Ce3. Therefore, this model is adopted in the present study:
"
C e3
6b ¼ 1:95 1 pbm
1=3 # ð33Þ
Cpp ¼
X 1 Cpw þ Cps 1þX 1þX
ð40Þ
– Density of the dispersed phase.
qp ¼ qsa ð1 þ XÞ ¼ qs ð1 dÞð1 þ XÞ
ð41Þ
where qsa is the apparent density of solids. The thermodynamics properties of water vapor and dry air are calculated as a function of temperature from formulas given by Reynolds [37].
where bm = 0.64 is the random close-packing solid volume fraction. 2.10. Boundary conditions 2.8. Particle–particle and particle–wall collisions In the present study the particles are traced one by one and the colliding partner is generated numerically based on the average solid properties in the cell enclosing the particle. The collision probability is calculated as proposed by Oesterlé and Petitjean [35]. If the collision takes place, the post-collision velocities are calculated using the hard sphere model described by Crowe et al. [32]. The particle trajectory after impact with the wall is greatly affected by the wall morphology. As a result of the surface roughness, the particle is generally hits a local surface slightly inclined to the flow direction by a small angle. The virtual-wall model developed by Sommerfeld [36] is adopted in the present study to calculate the roughness contribution angle. The particle–particle and particle–wall collision models used in the present study are presented in details in Ref. [26]. The particle–particle and particle–wall restitution coefficients are taken to be 0.75 and 0.9, respectively. The friction coefficient is
In the present study, there are four types of boundary conditions. At inlet, the gas velocity, temperature, mass fraction of water vapor and turbulent kinetic energy and its dissipation rate are specified. In compressible flow computations, the gas mass flow rate is specified at inlet instead of the gas velocity. The inlet velocity profile is assumed to be uniform and the turbulent kinetic energy and its dissipation rate are calculated by: 3=2
kin ¼ 0:04U 2in ;
ein ¼ C3=4 l
kin 0:03D
ð42Þ
At outlet the gradient of all variables in the flow direction is set to zero, except the axial gas velocity which is corrected to satisfy the overall mass balance. At the wall, the no-slip boundary conditions are imposed for the momentum equation, while for the energy equation; adiabatic wall boundary condition is considered. At the centerline, the symmetric boundary conditions are applied.
164
S.M. El-Behery et al. / Computers & Fluids 68 (2012) 159–167
400
Exp., Ref. [1] Modified Weber [7] Baeyens et al. [1] Modified Ranz-Marshall [11]
Tg (K)
380
Table 3 Coefficient of determination of the tested heat and mass transfer correlations.
Debrand [31] Frantz [9]
Coefficient of determination, R2
Correlation
De Brandt [1]
Frantz [27] De Brandt [1,11] Debrand [20] Baeyens et al. [1] Modified Ranz–Marshall [11] Modified Weber [17]
360
340
Gas temperature
Solid water content
0.4836 0.8172 0.7237 0.9927a 0.9795 0.9227
0.1268 0.7121 0.4163 0.9866a 0.8186 0.8916
320 a
300
0
5
10
15
20
The largest R2.
25
x (m)
(a) Gas temperature
400
0.4
Exp., Ref. [1]
390
Present 2D Num.
0.3
Tg (K)
X (kg water /kg solid )
380
1D Num., Ref. [11] 2D Num., Ref. [15]
370
3D Num., Ref. [18]
360
0.2 350 340
0.1
330 0.0
0
5
10
15
20
0
5
25
Baeyens et al. [1]
o
Rocha [43]
Case I
Case II
1.25 25.0 190 700 PVC 140 1195 980 400 0.26 0.125 12.81
1.25 25.0 190 700 PVC 180 1116 980 399 0.4 0.125 12.9
0.0525 4.0 100 150 Sand 380 2622 799.7 382.4 0.0468 0.0381 0.03947
20
25
1.58
1.85
0.00474
300 0.01
300 0.01
312.9 0.0469
2.11. Solution procedure Finite volume discretizations using the hybrid scheme for all variables, expect the density which is interpolated using the first order upwind scheme, are applied. The iterative solution based on the SIMPLE algorithm of Patankar [38] is used for the solution of the gas phase with an extended technique to compressible flow according to Karki [39]. The equations of motion of each particle along with its temperature equation are integrated using fourth
0.25
X (kg water /kg solid )
Table 2 Flow conditions and physical properties of drying test cases.
Solid mass flow rate, ms (kg/s) Solid inlet temperature, Tp,in (K) Mass fraction of water vapor at inlet Y H2 O;in
25
0.30
Fig. 2. Comparisons between present predictions using different correlations and experimental data of Baeyens et al. [1] (Case I).
o
20
(a) Gas temperature
(b) Solid water content
Gas mass flow rate, mg (kg/s)
15
x (m)
x (m)
Pipe diameter, D (m) Pipe length, L (m) Computational grid (radial axial) Particle material Particle diameter, dp (lm) Solid density, qs (kg/m3) Solid specific heat, Cps (J/kg/K) Inlet gas temperature, Tg,in (K) Inlet water content of solid, Xin Critical water content, Xcr
10
0.20
0.15
0.10
0.05
0
5
10
15
x (m)
(b) Solid water content Fig. 3. Comparisons between present predictions and predictions of Refs. [11,15,18] with experimental data of Ref. [1] (Case II).
order Runge–Kutta method. The grid used in the present work is non uniform in the radial direction. Thus, the grid is very fine near the pipe wall and gradually expanded to the pipe centerline. The computational grid is selected based on grid independence study and the dimensionless wall distance, y+ is less than unity in all the tested cases. Furthermore, the number of computational parcels and the time step used in the integration of Eqs. (6)–(10) are selected to insure independent results. The FORTRAN code used in the present study was originally developed by El-Behery et al. [40,41] to calculate gas–solid flow in curved ducts. The code is extended by El-Behery et al. [26] to include four-way coupling. Heat transfer and the compressibility effects were implemented by ElBehery et al. [42]. In the present study the code is modified to include mass transfer between both phases. The solution procedure for the fluid and particulate phases is as follows:
165
S.M. El-Behery et al. / Computers & Fluids 68 (2012) 159–167
390
322 Exp., Ref. [43] Present 2D Num.
380
320
1D Num., Ref. [11]
Tp (K)
Tg (K)
2D Num., Ref. [16] 3D Num., Ref. [18]
370
318
316
360 314
350
312 0
1
2
3
4
0
1
x (m)
2
3
4
x (m)
(a) Gas temperature
(b) Solid temperature
0.054
0.05
X (kg water /kg solid )
YH2O (kg water /kg air )
0.04
0.052
0.050
0.03
0.02
0.048 0.01
0.046
0.00 0
1
2
3
4
0
1
2
3
x (m)
x (m)
(c) Air water content
(d) Solid water content
4
Fig. 4. Comparisons between present prediction and predictions of Refs. [11,16,18] with experimental data of Ref. [43].
1. A converged solution of gas phase is calculated without source term of the dispersed phase and with gas void fraction, a of unity. The solution is converged when the normalized residuals are less than 0.001. 2. By numerically integrating the transitional and rotational equation of motion for each parcel, a large number of discrete parcels are traced through the flowfield. In the first iteration, the particle motion is obtained without particle–particle collision and information is stored for each cell to calculate the collision probability in the next iteration. 3. The void fraction for dispersed phase, b and for gas phase, a as well as the source terms are calculated. 4. The gas flowfield is recalculated taking into account the source terms and void fractions resulting in step 3. 5. Repeat steps 2–4 until the maximum change in the axial gas velocity between two successive coupled iterations is less than 0.001 of the mean gas velocity.
by Skuratovsky et al. [15,16], Mezhericher et al. [17] and Jamaleddine and Ray [18], the wall temperature is varied linearly form 325 K at inlet to 320 K at outlet. Other conditions for this test case are given in Table 2. For direct comparison with experimental data, mass weighted average technique is used to obtain the average value of various solution properties (except temperature), Refs. [15,16], as follows:
/¼
2p
RR
2p
qU/r dr qUr dr 0
R0 R
ð43Þ
The average temperature of each phase (the mean bulk temperature) can be calculated as given [15,16] as:
T¼
2p
RR
qUC p Tr dr 2p 0 qUC p r dr 0
RR
ð44Þ
3. Results and discussions The present model was thoroughly validated by El-Behery et al. [26,42] for hydrodynamics and thermal fields. The results presented herein concern only with heat and mass transfer processes. Six popular correlations for heat and mass transfer coefficients are tested in the present study (i.e., Eqs. (19)–(25)). Fig. 2 shows a comparison between present model predictions and experimental data of Baeyens et al. [1] (Case I) using different correlations for heat and mass transfer coefficients. Following the suggestion given
Fig. 2 indicates that Debrand, modified Weber and modified Ranz–Marshall correlations under-predict the gas temperature and the solid water content at the dryer outlet. On the other hand, Frantz and De Brandt correlation over-predict them. Baeyens et al. [1] correlation gives the best overall performance. For quantitative assessment of these correlations, the coefficient of determination, R2 is calculated for each correlation and presented in Table 3. The coefficient of determination, R2 is a measure of the goodness of a model and it can be defined as the ratio between residual sum of
166
S.M. El-Behery et al. / Computers & Fluids 68 (2012) 159–167
X (kg water /kg solid )
0.4 Present 2D, r = 0.0 m Present 2D, r = 0.62 m Exp., Ref. [1]
0.3
3D Ref. [17], r = 0.0 m 3D Ref. [17], r = 0.62 m
0.2
0.1
0.0
0
5
10
15
20
25
20
25
x (m)
(a) Solid water content 400
Tg (K)
380
360
340
320
0
10
5
15
Furthermore, the present predictions are in close agreement with the one-dimensional model. The main advantage of 2D and 3D models is that the cross-sectional distribution of flow parameters can be obtained. To assess the accuracy of the present model for the radial variation of flow properties, the present predictions are compared with the 3D Eulerian–Lagrangian predictions of Mezhericher et al. [17] using FLUENT CFD code, as shown in Fig. 5. The results are presented for two radial positions namely; pipe center (r = 0.0 m) and near periphery (r = 0.62 m). The experimental data of Baeyens et al. [1] are also presented in the same figure for direct comparison. This figure indicates that gas temperature is lower near the pipe wall due to the lower wall temperature. This results in low drying rate and hence higher solid moisture content near the pipe wall. The figure indicates also that the present results are more realistic than those of Mezhericher et al. [17] as it compared with experimental data. This may be due to the relatively coarse grid used by Mezhericher et al. [17]. Another improvement in the present model over that of Mezhericher et al. [17] is the use of variable turbulent Prandtl number. This parameter has a great effect on the heat transfer rate through the wall as reported by El-Behery et al. [26]. In order to understand the relative influence of turbulence dispersion and particle–particle collision, numerical simulation is carried out neglecting either of these parameters. A comparison between the complete model prediction and the predictions without either turbulence dispersion or particle–particle collision for Baeyens et al. [1] (Case I) is shown Fig. 6. It can be seen from this figure that the particle–particle collision is less important for this case. This can be attributed to the low solid void fraction
x (m)
(b) Gas temperature Fig. 5. Comparisons between present prediction and predictions of Ref. [17] with experimental data of Ref. [1] (Case I) at different radial positions.
R2 ¼ 1 Pi¼1 n
360
2
ð/i;num /i;exp Þ
i¼1 ð/i;exp
Complete model Without particle-particle collision Without turbulence dispersion
380
/exp Þ2
ð45Þ
The results presented in Table 3 indicate that the best prediction (the highest R2) is obtained when Baeyens et al. correlation is applied. In addition, modified Ranz–Marshall and modified Weber correlations can be applied to the cases involving heat transfer only. Overall, however, the correlation proposed by Baeyens et al. [1] is recommended in the present study for the heat and mass transfer coefficients. Figs. 3 and 4 present comparisons between present predictions and experimental data of Baeyens et al. [1] (Case II) and Rocha [43], respectively. The experimental data of Rocha are reported in Refs. [11,14–18]. The pipe wall temperature for Rocha test case is falling linearly from 360 K at inlet to 354 K at outlet. Other conditions for these test cases are given in Table 2. The figures also present direct comparison between the present predictions and other computations from literature. The selected numerical results for comparison are the one-dimensional prediction of Levy and Borde [11], two-dimensional numerical results of Skuratovsky et al. [15,16] and FLUENT three-dimensional predictions of Jamaleddine and Ray [18]. These models are based on Eulerian–Eulerian approach. The figures show that the present model agrees well with experimental data for both cases. In addition, the present model predicts the temperature and water content better than other two-dimensional and three-dimensional Eulerian–Eulerian models.
340
320
0
5
10
15
20
25
20
25
x (m)
(a) Gas temperature 0.4
X (kg water /kg solid )
Pn
Exp., Ref. [1]
Tg (K)
squares and explained sum of squares (i.e., total variation in experimental data), see Cameron and Windmeijer [44]. R2 can be calculated by:
400
0.3
0.2
0.1
0.0
0
5
10
15
x (m)
(b) Solid water content Fig. 6. Effect of particle–particle collision and turbulence dispersion on the model predictions.
S.M. El-Behery et al. / Computers & Fluids 68 (2012) 159–167
(<1.3 104) in this case. On the other hand, turbulence dispersion has a great effect on the gas temperature and solid water content. 4. Conclusion A steady-state two-dimensional four-way coupling Eulerian– Lagrangian model is presented. The model takes into account momentum, heat and mass exchange between both phases. Several correlations for heat transfer coefficient are tested in the present study. The model is validated with experimental data from previous investigators under different conditions. The present results are also compared to other published computations. It is found that the correlation proposed by Baeyens et al. [1] performs better than other correlations. The comparisons with experimental data showed that the present model is able to predict heat and mass transfer in gas flow with a good accuracy. In addition, the present model performs better than other models available in the literature. It is found also that neglecting the particle–particle collision or turbulence dispersion results in lower heat transfer and drying rates. References [1] Baeyens J, Gauwbergen DV, Vinckier I. Pneumatic drying: the use of large-scale data in a design procedure. Powder Technol 1995;83:139–48. [2] Borde I, Levy A. Pneumatic and flash drying. In: Mujumdar AS, editor. Handbook of industrial drying. New York: CRC Press; 2006. [3] Indrato A, Halim Y, Partoputro P. Pneumatic drying of solid particle: experimental and model comparison. Exp Heat Transfer 2007;20:277–87. [4] Williams-Gardner. Industrial drying. London, UK: Leonard Hill; 1971. [5] Moyers CG, Baldwin GW. Psychrometry, evaporative cooling, and solids drying. In: Perry RH, Green DW, Maloney JO, editors. Perry’s chemical engineers’ handbook. McGraw-Hill, Inc.; 1997. [6] Thorpe GR, Wint A, Coggan GC. The mathematical modelling of industrial pneumatic driers. Trans Inst Chem Eng 1973;51:339–48. [7] Kemp IC, Bahu RE, Pasley HS. Model development and experimental studies of vertical pneumatic conveying dryers. Dry Technol 1994;12:1323–40. [8] Kemp IC, Oakley DE. Simulation and scale-up of pneumatic conveying and cascading rotary dryers. Dry Technol 1997;15:1699–710. [9] Radford RD. A model of particulate drying in pneumatic conveying systems. Powder Technol 1997;93:109–26. [10] Kemp IC, Oakley DE, Bahu RE. Computational fluid dynamics modelling of vertical pneumatic conveying dryers. Powder Technol 1991;65:477–84. [11] Levy A, Borde I. Two-fluid model for pneumatic drying of particulate materials. Dry Technol 2001;19:1773–88. [12] Pelegrina AH, Crapiste GH. Modelling the pneumatic drying of food particles. J Food Eng 2001;48:301–10. [13] El-Behery SM, El-Askary WA, Ibrahim KA, Hamed MH. Porous particles drying in a vertical upward pneumatic conveying dryer. Int J Aerospace Mech Eng 2011;5:110–25. [14] Narimatsu CP, Ferreira MC, Freire JT. Drying of coarse particles in a vertical pneumatic conveyor. Dry Technol 2007;25:291–302. [15] Skuratovsky I, Levy A, Borde I. Two-fluid two-dimensional model for pneumatic drying. Dry Technol 2003;21:1645–68. [16] Skuratovsky I, Levy A, Borde I. Two-dimensional numerical simulations of the pneumatic drying in vertical pipes. Chem Eng Process 2005;44:187–92. [17] Mezhericher M, Levy A, Borde I. Three-dimensional modelling of pneumatic drying process. Powder Technol 2010;203:371–83. [18] Jamaleddine TJ, Ray MB. Numerical simulation of pneumatic dryer using computational fluid dynamics. Ind Eng Chem Res 2010;49:5900–10.
167
[19] Jamaleddine TJ, Ray MB. Drying of sludge in a pneumatic dryer using computational fluid dynamics. Dry Technol 2011;29:308–22. [20] NamKung W, Cho M. Pneumatic drying of iron ore particles in a vertical tube. Dry Technol 2004;22:877–91. [21] Kaensup W, Kulwong S, Wongwises S. A small-scale pneumatic conveying dryer of rough rice. Dry Technol 2006;24:105–13. [22] Kaensup W, Kulwong S, Wongwises S. Comparison of drying kinetics of paddy using a pneumatic dryer with and without a cyclone. Dry Technol 2006;24:1039–45. [23] Chagras V, Oesterlé B, Boulet P. On heat transfer in gas–solid pipe flows: effects of collision induced alterations of the flow dynamics. Int J Heat Mass Transfer 2005;48:1649–61. [24] Boulet P, Moissette S. Influence of the particle–turbulence modulation modeling in the simulation of a non-isothermal gas–solid flow. Int J Heat Mass Transfer 2002;45:4201–16. [25] Mansoori Z, Saffar-Avval M, Tabrizi HB, Ahmadi G. Modeling of heat transfer in turbulent gas–solid flow. Int J Heat Mass Transfer 2002;45: 1173–84. [26] El-Behery SM, El-Askary WA, Hamed MH, Ibrahim KA. Hydrodynamic and thermal fields analysis in gas–solid two-phase flow. Int J Heat Fluid Flow 2011;32:740–54. [27] Mansoori Z, Saffar-Avval M, Tabrizi HB, Dabir B, Ahmadi G. Inter-particle heat transfer in a riser of gas–solid turbulent flows. Powder Technol 2005;159:35–45. [28] Launder BE, Sharma BI. Application of the energy-dissipation model of turbulence to the calculation of flow near a spinning disc. Lett Heat Mass Transfer 1974;1:131–8. [29] Kays WM. Turbulent Prandtl number – where are we? J Heat Transfer 1994;116:284–95. [30] Shuen JS, Chen LD, Faeth GM. Evaluation of a stochastic model of particle dispersion in a turbulent round jet. AIChE J 1983;29:167–70. [31] Debrand S. Heat transfer during a flash drying process. Ind Eng Chem Process Des Develop 1974;13:396–404. [32] Crowe C, Sommerfeld M, Tsuji Y. Multiphase flow with droplets and particles. Florida, USA: CRC Press; 1998. [33] Lain S, Sommerfeld M. Turbulence modulation in dispersed two-phase flow laden with solids from a Lagrangian perspective. Int J Heat Fluid Flow 2003;24:616–25. [34] Zhang YH, Reese JM. Gas turbulence modulation in a two-fluid model for gas– solid flows. AICHE J 2003;49:3048–65. [35] Oesterlé B, Petitjean A. Simulation of particle-to-particle interactions in gas– solid flows. Int J Multiphase Flow 1993;19:199–211. [36] Sommerfeld M. Modelling of particle–wall collisions in confined gas–particle flows. Int J Multiphase Flow 1992;18:905–26. [37] Reynolds WC. Thermodynamic properties in SI: graphs, tables and computational equations for forty substances. Published by the Department of Mechanical Engineering. Stanford, CA 94305: Stanford University; 1979. [38] Patankar SV. Numerical heat transfer and fluid flow. New York, USA: McGrawHill; 1980. [39] Karki KC. A calculation procedure for viscous flows at all speeds in complex geometries. Ph.D. thesis. University of Minnesota; 1986. [40] El-Behery SM, Hamed MH, El-Kadi MA, Ibrahim KA. CFD prediction of air–solid flow in 180° curved duct. Powder Technol 2009;191:130–42. [41] El-Behery SM, Hamed MH, Ibrahim KA, El-Kadi MA. CFD evaluation of solid particles erosion in curved ducts. ASME J Fluids Eng 2010;132: 071303–13. [42] El-Behery SM, El-Askary WA, Hamed MH, Ibrahim KA. Numerical and experimental studies of heat transfer in particle-laden gas flows through a vertical riser. Int J Heat Fluid Flow 2012;33:118–30. [43] Rocha SCS. Contribution to the study of pneumatic drying: simulation and influence of the gas–particle heat transfer coefficient. Ph.D. thesis. Sao Paulo: EPUSP, Sao Paulo University; 1988. [44] Cameron AC, Windmeijer FAG. An R-squared measure of goodness of fit for some common nonlinear regression models. J Economet 1997;77:329–42.