**2. Numerical models**

86 Mass Transfer in Chemical Engineering Processes

modeling for the simulation of the drying mechanism in commercial dryers demonstrated the CFD capabilities and usefulness for the design and understanding of drying equipment. In a recent paper by Jamaleddine and Ray (2010)[3], the authors presented a comprehensive review on the application of CFD for the design, study, and evaluation of lab-scale and industrial dryers. The use of different numerical methods such as the finite element, finite volume, and finite difference were fully discussed. Numerical models such as the Eulerian-Eulerian and Eulerian-Lagrangian, used for gas–solid multiphase flow systems were also discussed along with their merits, disadvantages, and the scope of their applicability. The application of Kinetic theory approach for granular flow was also discussed. The authors pointed out some of the merits and shortcomings of CFD methods in general, and the drying application, in particular. They argued that a key advantage of CFD methods in evaluating drying systems is that it makes it possible to evaluate geometric changes (different feed point layouts such as multiple entry points) and operating conditions with much less time (faster turnaround time) and expense (flexibility to change design parameters without the expense of hardware changes) than would be involved in laboratory testing. A second advantage is that CFD provides far more detailed output information (suited for trouble-shooting) and far better understanding of the dryer performance than can be obtained in a laboratory environment. By interpreting graphical predictions from a CFD solution, local conditions of all phases in the drying chamber can be evaluated and crucial

information related to the dispersion of particulate material can be gathered.

agglomeration, and sintering are difficult to account for.

Despite the fact that CFD methods can offer valuable information and a great deal of insight of the process, the use of CFD methods requires considerable expertise. Lack of in-depth knowledge of the CFD methods and insufficient proficiency in utilizing commercial CFD software packages are major concerns for implementing CFD solutions in unknown and unconventional systems. In addition, CFD models have inherent limitations and challenges. Massah et al. (2000)[4] indicated some of the computational challenges of CFD modeling in the drying applications of granular material as follows. First of all, most processes involve solids with irregular shapes and size distribution, which might not be easily captured by some models. Second, Eulerian-Eulerian CFD methods rely on the kinetic theory approach to describe the constituent relations for solids viscosity and pressure, which are based on binary collisions of smooth *spherical* particles and do not account for deviations in shape or size distribution. Finally, very little is known about the turbulent interaction between different phases; thus, CFD models might not have the ability of presenting the associated drag models for a specific case study especially when solids concentration is high. In addition to the above, note that CFD simulations of three-dimensional geometries are computationally demanding and might be costly and although in some cases, the computational effort can be reduced by modeling a two-dimensional representation of the actual geometry (mostly for axisymmetric systems), the realistic behaviour of the simulated system might not be fully captured. Some geometrical systems cannot be modeled using the above simplification and thus, the computational effort becomes a must. This argument also applies to models adopting the Eulerian-Lagrangian formulation for *dense systems* which determine the trajectories of particles as they travel in the computational domain. In addition, formulas describing cohesion and frictional stresses within solids assembly are also not well established in these models. Finally, changes in particle size due to attrition, Multiphase flow models have improved substantially during the past years due to a better understanding of the physical phenomena occurring in multiphase flow systems. An extensive research has also led to a better understanding of the kinetic theory for granular flow and therefore, better implementation of the mathematical formulations pertaining to the flow, heat, and mass transfer mechanisms occurring in multiphase flow systems. The present numerical models for multiphase flows incorporate two approaches: the Eulerian-Eulerian approach, and the Eulerian-Lagrangian approach. A decision on whether the Eulerian-Eulerian or Eulerian-Lagrangian formulation of the governing equations is to be used should be made prior to the numerical solution, simply because each formulation has its limitations and constraints. Numerical predictions obtained from each formulation are not identical, and the choice of a convenient formulation for a specific model relies on whether a dense or dilute system is being considered and the objectives of the numerical study. For instance, if the objective of the numerical model is to follow the trajectories of individual particles, then the Eulerian-Lagrangian formulation appears more convenient for a dilute system (volume fraction of 1% and less). However, for a dense system, this approach is computationally expensive and time consuming and requires powerful and high-speed computers. On the contrary, the Eulerian-Eulerian formulation can handle both dense and dilute systems; however, it cannot predict the local behavior of particles in the flow field.

The theory behind the Eulerian-Eulerian approach is based on the macroscopic balance equations of mass, momentum, and energy for both phases. Eulerian models assume both phases as two interpenetrating continuum (Enwald et al., 1996)[5] and permit the solution of the Navier-Stokes equations with the assumption of incompressibility for both the gas and dispersed phases. The gas phase is the primary or continuous phase while the solid phase is termed as the dispersed phase. Both phases are represented by their volume fractions and are linked through the drag force in the momentum equation as given by Wen and Yu[6] correlation for a dilute system, Ergun[7] correlation for a dense system, and Gidaspow et al.[8], which is a combination of both correlations for transition and fluctuating systems. An averaging technique for the field variables such as the gas and solid velocities, solid volume

Numerical Simulation of Pneumatic and Cyclonic Dryers Using Computational Fluid Dynamics 89

using local flow parameters and gas properties, which is difficult to achieve using a continuum or steady-state model. The total number of particles is tractable from a computational point of view and modeling particle–particle and particle–wall interactions can be achieved with a great success. For additional information on the actual form of the conservation equations used in this approach, refer to Strang and Fix[15] and Gallagher[16]. In order to extend the applicability of single phase equations to multiphase flows, the volume fraction of each phase is implemented in the governing equations as was mentioned earlier. In addition, solids viscosities and stresses need to be addressed. The governing equations satisfying single phase flow will not be sufficient for flows where inter-particle interactions are present. These interactions can be in the form of collision between adjacent particles as in the case of a dilute system, or contact between adjacent particles in the case of dense systems. In the former, dispersed phase stresses and viscosities play a crucial role in the overall velocity and concentration distribution in the physical domain. The crucial factor attributed to this random distribution of particles in these systems is the gas phase turbulence. In cases where particles are light and small, turbulence eddies dominate the particles movement and the interstitial gas acts as a buffer that prevents collision between particles. However, in the case of heavy and large diameter particles (150 mm and higher), particle inertia is sufficient to carry them easily through the intervening gas film, and interactions occur by direct collision. Therefore, solids viscosities and stresses cannot be neglected, and the single phase fundamental equations need to be adjusted to account for

In the previous section, it was mentioned that each phase is represented by its volume fraction with respect to the total volume fraction of all phases present in the computational domain. For the sake of simplicity, let us develop these formulations for a binary system of two phases, a gas phase represented by *g*, and a solid phase represented by *s*. Accordingly, the mass conservation equation for each phase *q*, such that *q* can be a gas*= g* or solid= *s* is:

 

*g g g g g g g g <sup>g</sup> g g*

*U U U P <sup>g</sup> <sup>t</sup>*

*s s s s s s s s s s s s*

*U UU P <sup>g</sup> <sup>t</sup>*

*q pq qq qq*

1

. Similarly, the momentum balance equations for both

(1)

 

 

(2)

(3)

*n*

*p U M*

(defined later) represents the mass transfer from the *p*th phase to the *q*th phase.

*sg gs gs vm*

 

 

*MU F*

*sg gs gs vm*

*MU F*

the secondary phase interaction as shown in the next section.

*t* 

> 

> >

**2.2 Hydrodynamic model equations** 

When *q=g*, *p=s*, *Mpq M M sg gs*

where *Mpq* 

phases are:

fraction, and solid granular temperature is adopted. With this approach, the kinetic theory for granular flow (KTGF) is adopted to describe the interfacial forces between the considered phases and between each of the phases and the boundaries of the computational domain. The KTGF is based on the flow of nonuniform gases primarily presented by Chapman and Cowling (1970)[9]. The model was then further matured through the work of Jenkins and Savage (1983)[10], Lun et al. (1984)[18], Ding and Gidaspow (1990)[11], Gidaspow et al. (1992)[8], and Gidaspow (1994)[12].

On the other hand, Lagrangian models, or discrete particle models, are derived from Newton's law of motion for the dispersed phase. This approach facilitates the ability to compute the trajectory (path) and motion of individual particles. The interactions between the particles are described by either a potential force (soft particle dynamics)[13] or by collision dynamics (hard particle dynamics)[14]. In the Lagrangian approach, the fluid phase is treated separately by solving a set of time averaged Navier-Stokes equations, whereas the dispersed phase is solved by tracking a large number of particles, bubbles, or droplets in the calculated flow field. By computing the temporal development of a sufficiently large sample of particles, ensemble average quantities describing system performance can be evaluated. Furthermore, using the Lagrangian approach, the dispersed phase can exchange mass, momentum, and energy with the fluid phase through a source term added to the conservation equations. These equations also account for the changes in volume fraction of each phase. As each individual particle moves through the flow field, its trajectory, mass, and heat transfer calculations are obtained from a force balance along with an updated local conditions of the continuous phase by solving mass and energy balance equations. Thus, external forces acting on the solid particle such as aerodynamic, gravitational, buoyancy, and contact due to collisions among the particles and between the particles and the domain boundaries, can be calculated simultaneously with the particle motion using local parameters of gas and solids.

Although the form of the Eulerian momentum equation can be derived from its Lagrangian equivalent by averaging over the dispersed phase, each model has its advantages and disadvantages depending on the objective of the study and the type of system used. With this and the above definitions in mind, we now discuss the merits and shortcomings of each formulation.

### **2.1 Merits and shortcomings of each approach**

Some of the advantages and disadvantages of the Eulerian and Lagrangian formulations are discussed in this section. Examples of their use for actual physical systems are also provided to facilitate and enhance our understanding of the subject and to direct the reader to the appropriate formulation for the problem at hand.

For modeling spray dryers, coal and liquid fuel combustion, and particle-laden flows, the Lagrangian description of the governing equations is more suitable because these systems are considered dilute; that is, they are characterized by low concentration of particles with solid volume fractions on the order of 1% or less. It was previously mentioned that this characteristic of particles density allows the tracking of particles trajectories at different locations in the computational domain with less computational effort than the case for a dense system. Predicting the particles trajectories is the main distinctive advantage of the Lagrangian technique over the Eulerian formulation. This in turn provides the opportunity to evaluate interactions between particles, fluids, and boundaries at the microscopic level conservation equations used in this approach, refer to Strang and Fix[15] and Gallagher[16]. In order to extend the applicability of single phase equations to multiphase flows, the volume fraction of each phase is implemented in the governing equations as was mentioned earlier. In addition, solids viscosities and stresses need to be addressed. The governing equations satisfying single phase flow will not be sufficient for flows where inter-particle interactions are present. These interactions can be in the form of collision between adjacent particles as in the case of a dilute system, or contact between adjacent particles in the case of dense systems. In the former, dispersed phase stresses and viscosities play a crucial role in the overall velocity and concentration distribution in the physical domain. The crucial factor attributed to this random distribution of particles in these systems is the gas phase turbulence. In cases where particles are light and small, turbulence eddies dominate the particles movement and the interstitial gas acts as a buffer that prevents collision between particles. However, in the case of heavy and large diameter particles (150 mm and higher), particle inertia is sufficient to carry them easily through the intervening gas film, and interactions occur by direct collision. Therefore, solids viscosities and stresses cannot be neglected, and the single phase fundamental equations need to be adjusted to account for the secondary phase interaction as shown in the next section.

### **2.2 Hydrodynamic model equations**

88 Mass Transfer in Chemical Engineering Processes

fraction, and solid granular temperature is adopted. With this approach, the kinetic theory for granular flow (KTGF) is adopted to describe the interfacial forces between the considered phases and between each of the phases and the boundaries of the computational domain. The KTGF is based on the flow of nonuniform gases primarily presented by Chapman and Cowling (1970)[9]. The model was then further matured through the work of Jenkins and Savage (1983)[10], Lun et al. (1984)[18], Ding and Gidaspow (1990)[11], Gidaspow et

On the other hand, Lagrangian models, or discrete particle models, are derived from Newton's law of motion for the dispersed phase. This approach facilitates the ability to compute the trajectory (path) and motion of individual particles. The interactions between the particles are described by either a potential force (soft particle dynamics)[13] or by collision dynamics (hard particle dynamics)[14]. In the Lagrangian approach, the fluid phase is treated separately by solving a set of time averaged Navier-Stokes equations, whereas the dispersed phase is solved by tracking a large number of particles, bubbles, or droplets in the calculated flow field. By computing the temporal development of a sufficiently large sample of particles, ensemble average quantities describing system performance can be evaluated. Furthermore, using the Lagrangian approach, the dispersed phase can exchange mass, momentum, and energy with the fluid phase through a source term added to the conservation equations. These equations also account for the changes in volume fraction of each phase. As each individual particle moves through the flow field, its trajectory, mass, and heat transfer calculations are obtained from a force balance along with an updated local conditions of the continuous phase by solving mass and energy balance equations. Thus, external forces acting on the solid particle such as aerodynamic, gravitational, buoyancy, and contact due to collisions among the particles and between the particles and the domain boundaries, can be calculated simultaneously with the particle motion using local

Although the form of the Eulerian momentum equation can be derived from its Lagrangian equivalent by averaging over the dispersed phase, each model has its advantages and disadvantages depending on the objective of the study and the type of system used. With this and the above definitions in mind, we now discuss the merits and shortcomings of each

Some of the advantages and disadvantages of the Eulerian and Lagrangian formulations are discussed in this section. Examples of their use for actual physical systems are also provided to facilitate and enhance our understanding of the subject and to direct the reader to the

For modeling spray dryers, coal and liquid fuel combustion, and particle-laden flows, the Lagrangian description of the governing equations is more suitable because these systems are considered dilute; that is, they are characterized by low concentration of particles with solid volume fractions on the order of 1% or less. It was previously mentioned that this characteristic of particles density allows the tracking of particles trajectories at different locations in the computational domain with less computational effort than the case for a dense system. Predicting the particles trajectories is the main distinctive advantage of the Lagrangian technique over the Eulerian formulation. This in turn provides the opportunity to evaluate interactions between particles, fluids, and boundaries at the microscopic level

al. (1992)[8], and Gidaspow (1994)[12].

parameters of gas and solids.

**2.1 Merits and shortcomings of each approach** 

appropriate formulation for the problem at hand.

formulation.

In the previous section, it was mentioned that each phase is represented by its volume fraction with respect to the total volume fraction of all phases present in the computational domain. For the sake of simplicity, let us develop these formulations for a binary system of two phases, a gas phase represented by *g*, and a solid phase represented by *s*. Accordingly, the mass conservation equation for each phase *q*, such that *q* can be a gas*= g* or solid= *s* is:

$$\frac{\partial}{\partial t}(\rho\_q \alpha\_q) + \nabla \cdot \left(\rho\_q \alpha\_q \overrightarrow{\mathbf{U}}\_q\right) = \sum\_{p=1}^n \overset{\bullet}{M}\_{pq} \tag{1}$$

where *Mpq* (defined later) represents the mass transfer from the *p*th phase to the *q*th phase. When *q=g*, *p=s*, *Mpq M M sg gs* . Similarly, the momentum balance equations for both phases are:

$$\begin{aligned} \frac{\partial}{\partial t} \left( \rho\_{\mathcal{S}} \boldsymbol{\alpha}\_{\mathcal{S}} \overrightarrow{\boldsymbol{U}}\_{\mathcal{S}} \right) + \nabla \cdot \left( \rho\_{\mathcal{S}} \boldsymbol{\alpha}\_{\mathcal{S}} \overrightarrow{\boldsymbol{U}}\_{\mathcal{S}} \overrightarrow{\boldsymbol{U}}\_{\mathcal{S}} \right) &= -\boldsymbol{\alpha}\_{\mathcal{S}} \nabla P + \nabla \cdot \overset{\bullet}{\boldsymbol{\tau}}\_{\mathcal{S}} + \rho\_{\mathcal{S}} \boldsymbol{\alpha}\_{\mathcal{S}} \overset{\bullet}{\mathcal{S}} + \\\\ \overrightarrow{\boldsymbol{\beta}}\_{\mathcal{S}^{\rm s}} + \boldsymbol{M}\_{\text{sg}} \overrightarrow{\boldsymbol{U}}\_{\mathcal{S}^{\rm s}} + \overrightarrow{\boldsymbol{F}}\_{\text{vm}} \end{aligned} \tag{2}$$
 
$$\begin{aligned} \frac{\partial}{\partial t} \left( \rho\_{s} \boldsymbol{\alpha}\_{s} \overrightarrow{\boldsymbol{U}}\_{s} \right) + \nabla \cdot \left( \rho\_{s} \boldsymbol{\alpha}\_{s} \overrightarrow{\boldsymbol{U}}\_{s} \overrightarrow{\boldsymbol{U}}\_{s} \right) &= -\boldsymbol{\alpha}\_{s} \nabla P\_{s} + \nabla \cdot \boldsymbol{\tau}\_{s} + \rho\_{s} \boldsymbol{\alpha}\_{s} \overset{\bullet}{\mathcal{S}} - \tag{3}$$

$$
\begin{aligned}
\vec{\beta}\_{\text{test}} &= \vec{\beta}\_{\text{test}} + \vec{\beta}\_{\text{test}} + \vec{\beta}\_{\text{test}} + \vec{\beta}\_{\text{test}} + \vec{\beta}\_{\text{test}} + \vec{\beta}\_{\text{test}} \\ &= \vec{\beta}\_{\text{ges}} - \vec{M}\_{\text{sig}} \vec{U}\_{\text{ges}} + \vec{F}\_{\text{sun}}
\end{aligned}
\tag{3}
$$

Numerical Simulation of Pneumatic and Cyclonic Dryers Using Computational Fluid Dynamics 91

system adopting the Eulerian formulation such that *q= g* for gas and *s* for solid, the volume fraction balance equation representing both phases in the computational domain can then be

1

(8)

 

1

In the case of collision between the particles in the solid phase, the kinetic theory for granular flow based on the work of Gidaspow et al. (1992)[8] dictates that the solid shear

> <sup>1</sup> <sup>2</sup> <sup>10</sup> 4 4 <sup>2</sup> 11 1

where *ss e* is a value between *0* and *1* dictating whether the collision between two solid particles is inelastic or perfectly elastic. When two particles collide, and depending on the material property, initial particle velocity, etc, deformation in the particle shape might occur. The resistance of granular particles to compression and expansion is called the solid bulk

*<sup>d</sup> g e d eg*

1 21 *ss s o*

*<sup>s</sup>* is the granular temperature which measures the kinetic energy fluctuation in the

: 3

 

*s s s ss s sg s*

 

*U*

where the first term on the right hand side (RHS) is the generation of energy by the solid stress tensor; the second term represents the diffusion of energy; the third term represents the collisional dissipation of energy between the particles; and the fourth term represents

2 3 *<sup>s</sup>* 

<sup>4</sup> <sup>2</sup> <sup>1</sup>

*<sup>s</sup> b s s s o ss dg e*

 

1

(10)

*<sup>g</sup>* (11)

*c* (12)

(13)

 

(9)

*ss s <sup>s</sup> <sup>s</sup> s o ss s s s ss o*

*<sup>s</sup>* can be represented by Equation (9) as follows:

*s ss o*

 

*e g*

96 1 5 5

*<sup>b</sup>* . According to Lun et al. (1984)[18] correlation, it is given by:

3

 

This parameter can be governed by the following conservation equation:

*s s ss s s s*

 

the energy exchange (transfer of kinetic energy) between the gas and solid phases.

*PI U k*

solid phase written in terms of the particle fluctuating velocity *c* as:

3 2

*t*

 

In addition, the solid pressure *Ps* is given by Gidaspow and Huilin (1998)[19] as:

*P e s s ss*

*n q q* 

given as:

where *<sup>q</sup> q V V*

viscosity

viscosity

where  such that *Ugs* is the relative velocity between the phases given by *U UU gs <sup>g</sup> <sup>s</sup>* . In the above equations, *gs* represents the drag force between the phases and is a function of the interphase momentum coefficient *Kgs* , the number of particles in a computational cell *Nd*, and the drag coefficient*CD* such that:

$$\begin{split} \overline{\mathcal{O}}\_{\mathcal{S}^{s}} &= \mathcal{K}\_{\mathcal{S}^{s}} \left( \overline{\mathcal{U}}\_{\mathcal{S}} - \overline{\mathcal{U}}\_{s} \right) \\ &= \mathcal{N}\_{d} \mathcal{C}\_{D} \frac{1}{2} \mathcal{P}\_{\mathcal{S}} \left( \mathcal{U}\_{\mathcal{S}} - \mathcal{U}\_{s} \right) \left| \overline{\mathcal{U}}\_{\mathcal{S}} - \overline{\mathcal{U}}\_{s} \right| A\_{surface} \\ &= \frac{6 \alpha\_{s}}{\pi d\_{s}^{3}} \mathcal{C}\_{D} \frac{1}{2} \mathcal{P}\_{\mathcal{S}} \left( \mathcal{U}\_{\mathcal{S}} - \mathcal{U}\_{s} \right) \left| \overline{\mathcal{U}}\_{\mathcal{S}} - \overline{\mathcal{U}}\_{s} \right| \frac{\pi d\_{s}^{2}}{4} \\ &= \frac{3}{4} \mathcal{C}\_{D} \frac{\alpha\_{s} \mathcal{P}\_{\mathcal{S}}}{d\_{s}} \Big( \mathcal{U}\_{\mathcal{S}} - \mathcal{U}\_{s} \Big) \left| \overline{\mathcal{U}}\_{\mathcal{S}} - \overline{\mathcal{U}}\_{s} \right| \end{split} \tag{4}$$

The form of the drag coefficient in Equation (4) can be derived based on the nature of the flow field inside the computational domain. Several correlations have been derived in the literature. A well established correlation that takes into consideration changes in the flow characteristics for multiphase systems is Ossen drag model presented in Skuratovsky et al. (2003)[17] as follows:

$$\begin{aligned} \mathbf{C}\_{d} &= \frac{64}{\pi \operatorname{Re}\_{s}} \Big( 1 + \frac{64}{2\pi} \Big) & \text{for } \operatorname{Re}\_{s} < 0.01\\ \mathbf{C}\_{d} &= \frac{64}{\pi \operatorname{Re}\_{s}} \Big( 1 + 10^{x} \Big) & \text{for } 0.01 < \operatorname{Re}\_{s} < 1.5\\ \mathbf{x} &= -0.883 + 0.906 \ln \left( \operatorname{Re}\_{s} \right) - 0.025 \ln^{2} \left( \operatorname{Re}\_{s} \right) & \text{(5)}\\ \mathbf{C}\_{d} &= \frac{64}{\pi \operatorname{Re}\_{s}} \Big( 1 + 0.138 \operatorname{Re}\_{s}^{0.992} \Big) \Big/ \operatorname{for} \, 1.5 < \operatorname{Re}\_{s} < 133\\ \ln \operatorname{C}\_{d} &= 2.0351 - 1.66 \ln \operatorname{Re}\_{s} + \ln^{2} \operatorname{Re}\_{s} - 0.0306 \ln^{3} \operatorname{Re}\_{s} \\ & \qquad \qquad \qquad \qquad \qquad \qquad \qquad \qquad \qquad \qquad \qquad \qquad \qquad \qquad \qquad \qquad \qquad \qquad \qquad \qquad \qquad \qquad \qquad \end{aligned} \tag{5}$$

The form of Reynolds number defined in Equation (5) is a function of the gas properties, the relative velocity between the phases, and the solid phase diameter. It is given by:

$$\text{Re}\_s = \frac{\rho\_\mathcal{g} \left| \overrightarrow{\mathcal{U}}\_\mathcal{S} - \overrightarrow{\mathcal{U}}\_s \right| d\_s}{\mu\_\mathcal{g}} \tag{6}$$

The virtual-mass force *F vm* in Equations (2) & (3) accounts for the force needed to accelerate the fluid surrounding the solid particle. It is given by:

$$
\overrightarrow{F}\_{vm} = \alpha\_s \rho\_g c\_{vm} (\frac{d\overrightarrow{U}\_g}{dt} - \frac{d\overrightarrow{U}\_s}{dt}) \tag{7}
$$

#### **2.3 Complimentary equations – granular kinetic theory equations**

When the number of unknowns exceeds the number of formulated equations for a specific case study, complimentary equations are needed for a solution to be possible. For a binary system adopting the Eulerian formulation such that *q= g* for gas and *s* for solid, the volume fraction balance equation representing both phases in the computational domain can then be given as:

$$\sum\_{q=1}^{n} \alpha\_q = 1\tag{8}$$

where *<sup>q</sup> q V V* 

90 Mass Transfer in Chemical Engineering Processes

is the relative velocity between the phases given by *U UU gs <sup>g</sup> <sup>s</sup>* .

of the interphase momentum coefficient *Kgs* , the number of particles in a computational cell

*NC U U U U A*

*s s D gg s g s*

*<sup>d</sup> C U UU U*

2 4

*s s*

2

2 3

in Equations (2) & (3) accounts for the force needed to accelerate

(7)

*g s d D g g s surface*

*C U UU U*

*D gs g s*

The form of the drag coefficient in Equation (4) can be derived based on the nature of the flow field inside the computational domain. Several correlations have been derived in the literature. A well established correlation that takes into consideration changes in the flow characteristics for multiphase systems is Ossen drag model presented in Skuratovsky et al.

*KU U*

*s g*

*s*

*d*

 

*C for*

*C for*

*C for*

 

Re 2

*s*

*s*

Re

the fluid surrounding the solid particle. It is given by:

*x*

*C*

Re

*s*

*for*

*d ss*

*d s*

*x d s*

40 Re 1000

 

*s*

relative velocity between the phases, and the solid phase diameter. It is given by:

*s*

 

**2.3 Complimentary equations – granular kinetic theory equations** 

0.792

ln 2.0351 1.66ln Re ln Re 0.0306ln Re

<sup>64</sup> 1 0.138Re 1.5 Re 133

*d ss s*

The form of Reynolds number defined in Equation (5) is a function of the gas properties, the

Re *<sup>g</sup> g s <sup>s</sup>*

( ) *<sup>g</sup> <sup>s</sup> vm s g vm dU dU F c dt dt*

When the number of unknowns exceeds the number of formulated equations for a specific case study, complimentary equations are needed for a solution to be possible. For a binary

 

*g*

 *U Ud* 

0.883 0.906ln Re 0.025ln Re

<sup>64</sup> 1 10 0.01 Re 1.5

64 64 1 Re 0.01

*gs g s gs*

3

*s*

*d*

3 4

represents the drag force between the phases and is a function

2

(4)

(5)

(6)

 

such that *Ugs*

In the above equations, *gs*

(2003)[17] as follows:

The virtual-mass force *F vm*

*Nd*, and the drag coefficient*CD* such that:

In the case of collision between the particles in the solid phase, the kinetic theory for granular flow based on the work of Gidaspow et al. (1992)[8] dictates that the solid shear viscosity *<sup>s</sup>* can be represented by Equation (9) as follows:

$$\mu\_s = \frac{10\,\rho\_s d\_s \sqrt{\theta\_s \pi}}{96a\_s \left(1 + e\_{ss}\right) g\_o} \left[1 + \frac{4}{5} a\_s g\_o \left(1 + e\_{ss}\right)\right]^2 + \frac{4}{5} a\_s \rho\_s d\_s \left(1 + e\_{ss}\right) g\_o \left(\frac{\theta\_s}{\pi}\right)^{\frac{1}{2}}\tag{9}$$

where *ss e* is a value between *0* and *1* dictating whether the collision between two solid particles is inelastic or perfectly elastic. When two particles collide, and depending on the material property, initial particle velocity, etc, deformation in the particle shape might occur. The resistance of granular particles to compression and expansion is called the solid bulk viscosity *<sup>b</sup>* . According to Lun et al. (1984)[18] correlation, it is given by:

$$
\mu\_{\rm b} = \frac{4}{3} \alpha\_s \rho\_s d\_s \mathbf{g}\_o \left( \mathbf{1} + \mathbf{e}\_{ss} \right) \left( \frac{\theta\_s}{\pi} \right)^{\frac{1}{2}} \tag{10}
$$

In addition, the solid pressure *Ps* is given by Gidaspow and Huilin (1998)[19] as:

$$P\_s = \alpha\_s \rho\_s \theta\_s \left[ 1 + 2 \left( 1 + e\_{ss} \right) \alpha\_s g\_o \right] \tag{11}$$

where *<sup>s</sup>* is the granular temperature which measures the kinetic energy fluctuation in the solid phase written in terms of the particle fluctuating velocity *c* as:

$$
\theta\_s = \mathfrak{c}^2 \bigvee \mathfrak{3} \tag{12}
$$

This parameter can be governed by the following conservation equation:

$$\begin{aligned} &\frac{3}{2} \left[ \frac{\partial}{\partial t} (\rho\_s a\_s \theta\_s) + \nabla \cdot \left( \rho\_s a\_s \overrightarrow{\mathbf{U}}\_s \theta\_s \right) \right] \\ &= \left( -P\_s \overrightarrow{I} + \tau\_s \right) \colon \nabla \cdot \overrightarrow{\mathbf{U}}\_s + \nabla \cdot \left( k\_s \nabla \theta\_s \right) - \gamma\_s - \mathfrak{D} \rho\_{\text{sg}} \theta\_s \end{aligned} \tag{13}$$

where the first term on the right hand side (RHS) is the generation of energy by the solid stress tensor; the second term represents the diffusion of energy; the third term represents the collisional dissipation of energy between the particles; and the fourth term represents the energy exchange (transfer of kinetic energy) between the gas and solid phases.

Numerical Simulation of Pneumatic and Cyclonic Dryers Using Computational Fluid Dynamics 93

Pr *p g cond*

The governing equation for the CDP is expressed in Equation (24). This equation can be used regardless of the method adopted to determine the critical moisture content. In cases when the critical moisture content is known, the FRP can then be expressed as shown in Equation (25) such that *X XX eq cr* . When the critical value is not known, Equation (26) can then be used as shown below. This equation was derived based on Fick's diffusion

( ) *c sat s CDR H O*

equation (convection-diffusion equation) is used as shown in Equation (27).

 

*t* 

*kM P T <sup>P</sup> M X*

2

 

(25)

(26)

(27)

*s s g*

*d RT RT* 

2 2 *v s Diffusion eq <sup>D</sup> <sup>M</sup> X X R* 

*eq FDR CDR cr eq*

*X X M M X X*

In order to obtain the water vapor distribution in the gas phase, the species transport

*<sup>g</sup> sg ggg gg g gg v g Y U Y DY M*

During the drying process, liquid water is removed and the particle density gradually

increases. With the assumption of no shrinkage, the particle density is expressed by:

*c k* 

*<sup>g</sup>*

*Sc*

*g v*

*D* 

(23)

(24)

*s cond s*

The diffusion coefficient *Dv* defined in the above equations is assumed to be constant. As the wet feed comes in contact with the hot carrier fluid, heat exchange between the phases occurs. In this stage, mass transfer is considered negligible. When the particle temperature exceeds the vaporization temperature, water vapor evaporates from the surface of the particle. This process is usually short and is governed by convective heat and mass transfer. This initial stage of drying is known as the constant or unhindered drying period (CDP). As drying proceeds, internal moisture within the particle diffuses to the surface to compensate for the moisture loss at that region, and diffusion mass transfer starts to occur. This stage dictates the transfer from the CDP to the second or falling rate drying period (FRP) and is designated by the critical moisture content. This system specific value is crucial in depicting which drying mechanism occurs; thus, it has to be accurate. However, it is not readily available and should be determined from experimental observations for different materials. An alternative approach that bypasses the critical value yet distinguishes the two drying periods is by drawing a comparison to the two drying rates. If the calculated value of diffusive mass transfer is greater than the convective mass transfer, then resistance is said to occur on the external surface of the particle and the CDP dominates. However, if the diffusive mass transfer is lower than the convective counterpart, then resistance occurs in

*Nu k <sup>h</sup> <sup>d</sup>*

the core of the particle and diffusion mass transfer dominates.

equation[22] for a spherical particle averaged over an elementary volume.

where

The diffusion coefficient for the solid phase energy fluctuation given by Gidaspow et al. (1992)[8] is:

$$k\_s = \frac{150\,\rho\_s d\_s \sqrt{\theta\_s \pi}}{384\left(1 + e\_{ss}\right) \text{g}\_o} \left[1 + \frac{6}{5} \alpha\_s \text{g}\_o \left(1 + e\_{ss}\right)\right]^2 + 2\alpha\_s^2 \rho\_s d\_s \text{g}\_o \left(1 + e\_{ss}\right) \left(\frac{\theta\_s}{\pi}\right)^{\frac{1}{2}}\tag{14}$$

The dissipation of energy fluctuation due to particle collision given by Gidaspow et al. (1992)[8] is:

$$\mathcal{I}\mathcal{Y}\_s = 3\alpha\_s^2 \rho\_s \mathcal{g}\_o \left(1 - e\_{ss}^2\right) \theta\_s \left[\frac{4}{d\_s} \left(\frac{\theta\_s}{\pi}\right)^{\frac{1}{2}} - \nabla \cdot \overrightarrow{\mathcal{U}}\_s\right] \tag{15}$$

The radial distribution function *<sup>o</sup> g* based on Ding and Gidaspow (1990)[11] model is a measure of the probability of particles to collide. For dilute phases, 1 *<sup>o</sup> g* ; for dense phases, *<sup>o</sup> g* .

$$\mathbf{g}\_o = \frac{3}{5} \left[ 1 - \left( \frac{\alpha\_s}{\alpha\_{s,\text{max}}} \right)^{\frac{1}{3}} \right]^{-1} \tag{16}$$

#### **2.4 Drying model equations – heat and mass transfer**

The conservation equation of energy (*q = g, s*) is given by:

$$\frac{\partial}{\partial t}(\rho\_q \alpha\_q H\_q) + \nabla \cdot \left(\rho\_q \alpha\_q \overrightarrow{\mathcal{U}}\_q H\_q\right) = -\alpha\_q \nabla P + \overleftarrow{\tau}\_q : \nabla \overrightarrow{\mathcal{U}}\_q + Q\_{pq} + \overset{\bullet}{\mathcal{M}}\_{pq} H\_q \tag{17}$$

By introducing the number density of the dispersed phase (solid in this case), the intensity of heat exchange between the phases is:

$$\mathcal{Q}\_{sg} = \mathcal{N}\_d \pi d\_s^2 h \left( T\_g - T\_s \right) = \frac{6\alpha\_s}{d\_s} h \left( T\_g - T\_s \right) = \frac{6\alpha\_s}{d\_s} m\_s c\_p \frac{dT\_s}{dt} \tag{18}$$

Many empirical correlations are available in the literature for the value of the heat- and mass-transfer coefficients. The mostly suitable for pneumatic and cyclone dryers are those given by Baeyens et al. (1995)[20] and De Brandt (1974)[21]. The Chilton and Colburn analogy for heat and mass-transfer are used as follows:

$$\mathrm{Nu}\_s = 0.15 \,\mathrm{Re}\_s \tag{19}$$

$$\mathrm{Nu}\_s = 0.16 \,\mathrm{Re}\_s^{1.3} \,\mathrm{Pr}^{0.67} \tag{20}$$

$$\text{Sh} = 0.15 \,\text{Re}\_s \tag{21}$$

$$\text{lSh} = 0.16 \,\text{Re}\_s^{1.3} \,\text{Sc}^{0.67} \tag{22}$$

where

92 Mass Transfer in Chemical Engineering Processes

The diffusion coefficient for the solid phase energy fluctuation given by Gidaspow et al.

The dissipation of energy fluctuation due to particle collision given by Gidaspow et al.

<sup>2</sup> 2 2 <sup>4</sup> 3 1 *<sup>s</sup> <sup>s</sup> s s s o ss s*

 

The radial distribution function *<sup>o</sup> g* based on Ding and Gidaspow (1990)[11] model is a measure of the probability of particles to collide. For dilute phases, 1 *<sup>o</sup> g* ; for dense phases,

*s g e U d*

,max

*s*

*q q* : *q pq qq q qq q q H UH P U Q M H pq <sup>q</sup> <sup>t</sup>*

By introducing the number density of the dispersed phase (solid in this case), the intensity

<sup>2</sup> 6 6 *s ss sg d s g s g s s p*

Many empirical correlations are available in the literature for the value of the heat- and mass-transfer coefficients. The mostly suitable for pneumatic and cyclone dryers are those given by Baeyens et al. (1995)[20] and De Brandt (1974)[21]. The Chilton and Colburn analogy

*dT Q N d h T T h T T mc <sup>d</sup> d dt* 

 

*ss s <sup>s</sup> <sup>s</sup> s o ss s s s o ss*

*<sup>d</sup> <sup>k</sup> g e dg e*

3 1 5 *<sup>s</sup> <sup>o</sup>*

*g*

384 1 5

**2.4 Drying model equations – heat and mass transfer**  The conservation equation of energy (*q = g, s*) is given by:

 

of heat exchange between the phases is:

for heat and mass-transfer are used as follows:

 

 

*e g*

*ss o*

<sup>2</sup> <sup>2</sup> <sup>150</sup> <sup>6</sup> 11 2 1

 

1

<sup>1</sup> <sup>1</sup> 3

> 

*s s*

(17)

0.15Re *Nus s* (19)

0.15Re *Sh <sup>s</sup>* (21)

1.3 0.67 0.16Re Pr *Nus s* (20)

1.3 0.67 0.16Re *Sh <sup>s</sup> Sc* (22)

 

(14)

<sup>1</sup> <sup>2</sup>

(15)

(16)

(18)

(1992)[8] is:

(1992)[8] is:

*<sup>o</sup> g* .

$$\mathbf{h} = \frac{\mathbf{N}\mu\_s k\_{cond}}{d\_s} \qquad \mathbf{Pr} = \frac{\mathbf{c}\_p \mu\_g}{k\_{cond}} \qquad \mathbf{Sc} = \frac{\mu\_g}{\rho\_g D\_v} \tag{23}$$

The diffusion coefficient *Dv* defined in the above equations is assumed to be constant. As the wet feed comes in contact with the hot carrier fluid, heat exchange between the phases occurs. In this stage, mass transfer is considered negligible. When the particle temperature exceeds the vaporization temperature, water vapor evaporates from the surface of the particle. This process is usually short and is governed by convective heat and mass transfer. This initial stage of drying is known as the constant or unhindered drying period (CDP). As drying proceeds, internal moisture within the particle diffuses to the surface to compensate for the moisture loss at that region, and diffusion mass transfer starts to occur. This stage dictates the transfer from the CDP to the second or falling rate drying period (FRP) and is designated by the critical moisture content. This system specific value is crucial in depicting which drying mechanism occurs; thus, it has to be accurate. However, it is not readily available and should be determined from experimental observations for different materials. An alternative approach that bypasses the critical value yet distinguishes the two drying periods is by drawing a comparison to the two drying rates. If the calculated value of diffusive mass transfer is greater than the convective mass transfer, then resistance is said to occur on the external surface of the particle and the CDP dominates. However, if the diffusive mass transfer is lower than the convective counterpart, then resistance occurs in the core of the particle and diffusion mass transfer dominates.

The governing equation for the CDP is expressed in Equation (24). This equation can be used regardless of the method adopted to determine the critical moisture content. In cases when the critical moisture content is known, the FRP can then be expressed as shown in Equation (25) such that *X XX eq cr* . When the critical value is not known, Equation (26) can then be used as shown below. This equation was derived based on Fick's diffusion equation[22] for a spherical particle averaged over an elementary volume.

$$\overset{\bullet}{M}\_{\rm CDR} = \frac{k\_c M}{d\_s} \left( \frac{P\_{\rm sat}(T\_s)}{RT\_s} - X\_{H2O} \frac{P}{RT\_{\rm g}} \right) \tag{24}$$

$$
\stackrel{\bullet}{M}\_{FDR} = \frac{X - X\_{eq}}{X\_{cr} - X\_{eq}} \stackrel{\bullet}{M}\_{CDR} \tag{25}
$$

$$\stackrel{\bullet}{M}\_{D\overline{\text{diffusion}}} = \frac{D\_v \rho\_s \pi^2}{R^2} \left(\overline{X} - X\_{eq}\right) \tag{26}$$

In order to obtain the water vapor distribution in the gas phase, the species transport equation (convection-diffusion equation) is used as shown in Equation (27).

$$\frac{\partial}{\partial t}(\rho\_{\mathcal{S}}\boldsymbol{a}\_{\mathcal{S}}\boldsymbol{Y}\_{\mathcal{S}}) + \nabla \cdot \left(\rho\_{\mathcal{S}}\boldsymbol{a}\_{\mathcal{S}}\overrightarrow{\boldsymbol{U}}\_{\mathcal{S}}\boldsymbol{Y}\_{\mathcal{S}}\right) = \nabla \cdot \left(\boldsymbol{a}\_{\mathcal{S}}\rho\_{\mathcal{S}}\boldsymbol{D}\_{v}\nabla\boldsymbol{Y}\_{\mathcal{S}}\right) + \overset{\bullet}{M}\_{\text{sg}}\tag{27}$$

During the drying process, liquid water is removed and the particle density gradually increases. With the assumption of no shrinkage, the particle density is expressed by:

Numerical Simulation of Pneumatic and Cyclonic Dryers Using Computational Fluid Dynamics 95

fluctuations in the volume fraction. When multiplied by the interchange coefficient *Kgs* , it

*<sup>g</sup> <sup>s</sup> dr <sup>s</sup> <sup>g</sup> gs s gs g D D*

such that *D DD <sup>g</sup> s ts*, *<sup>g</sup>* for Tchen Theory of multiphase flow (FLUENT 6.3 User's guide)[23]. The generation of turbulence kinetic energy due to the mean velocity gradients *Gk*,*<sup>g</sup>* is

, , :

,

*tg g*

 *C* 

The Reynolds stress tensor defined in Equation (13) for the continuous phase is based on the

, , <sup>2</sup>

 *g g g g g g g tg k U I UU g g tg <sup>g</sup> <sup>g</sup>* 

Time and length scales that characterize the motion of solids are used to evaluate the dispersion coefficients, the correlation functions, and the turbulent kinetic energy of the particulate phase. The characteristic particle relaxation time connected with inertial effects

> 1 , *<sup>s</sup> F sg g <sup>s</sup> gs V*

*g K C* 

 , , <sup>2</sup> <sup>1</sup> *t g*

> , , *sg t g t g*

*C*

*U L* 

 

*G UU U kg tg*

 *gg g* 

2

*g*

*g k*

> 

 

serves as a correction to the momentum exchange term for turbulent flows:

*U*

is defined in Equation (33). This velocity results from turbulent

*T*

*<sup>t</sup>*,*<sup>g</sup>* given in the above equation is written in terms of the turbulent

(33)

(34)

(35)

(36)

*T*

(37)

(38)

(39)

The drift velocity *Udr*

computed from:

The turbulent viscosity

kinetic energy of the gas phase as:

Boussinesq hypothesis[24] given by:

where

3

acting on a particulate phase is defined as:

 

**2.5.2 Dispersed phase turbulence equations** 

 

 

*t sg*

The Lagrangian integral timescale calculated along particle trajectories is defined as:

$$\rho\_s = \frac{\rho\_{H\_2O(l)} \rho\_{ds}}{X \left(\rho\_{ds} - \rho\_{H\_2O(l)}\right) + \rho\_{H\_2O(l)}} \tag{28}$$

### **2.5 Turbulence model equations**

To describe the effects of turbulent fluctuations of velocities and scalar quantities in each phase, the *k* multiphase turbulent model can be used for simpler geometries. Advanced turbulence models should be used for cases with swirl and vortex shedding (*RANS*, *k* ). In the context of gas-solid models, three approaches can be applied (FLUENT 6.3 User's guide)[23]: (1) modeling turbulent quantities with the assumption that both phases form a mixture of density ratio close to unity (mixture turbulence model); (2) modeling the effect of the dispersed phase turbulence on the gas phase and vice versa (dispersed turbulence model); or (3) modeling the turbulent quantities in each phase independent of each other (turbulence model for each phase). In many industrial applications, the density of the solid particles is usually larger than that of the fluid surrounding it. Furthermore, modeling the turbulent quantities in each phase is not only complex, but also computationally expensive when large number of particles is present. A more desirable option would then be to model the turbulent effect of each phase on the other by incorporating *source* terms into the conservation equations. This model is highly applicable when there is one primary phase (the gas phase) and the rest are dispersed dilute secondary phases such that the influence of the primary phase turbulence is the dominant factor in the random motion of the secondary phase.

### **2.5.1 Continuous phase turbulence equations**

In the case of multiphase flows, the standard *k* model equations are modified to account for the effect of dispersed phase turbulence on the continuous phase as shown below:

$$\frac{\partial}{\partial t} \left( \rho\_{\mathcal{S}} \boldsymbol{\alpha}\_{\mathcal{S}} k \right) + \nabla \cdot \left( \rho\_{\mathcal{S}} \boldsymbol{\alpha}\_{\mathcal{S}} \overrightarrow{\boldsymbol{L}}\_{\mathcal{S}} k \right) = \nabla \cdot \left( \boldsymbol{\alpha}\_{\mathcal{S}} \frac{\mu\_{\boldsymbol{t}, \mathcal{S}}}{\sigma\_{\boldsymbol{k}}} \nabla k\_{\mathcal{S}} \right) + \boldsymbol{\alpha}\_{\mathcal{S}} \mathbf{G}\_{\boldsymbol{k}, \mathcal{S}} - \boldsymbol{\alpha}\_{\mathcal{S}} \rho\_{\mathcal{S}} \boldsymbol{\varepsilon}\_{\mathcal{S}} + \boldsymbol{\alpha}\_{\mathcal{S}} \rho\_{\mathcal{S}} \Pi\_{\boldsymbol{k}\_{\mathcal{S}}} \tag{29}$$

and

$$\begin{split} \frac{\partial}{\partial t} \Big( \rho\_{\mathcal{S}} a\_{\mathcal{S}} \boldsymbol{\varepsilon}\_{\mathcal{S}} \Big) + \nabla \cdot \Big( \rho\_{\mathcal{S}} a\_{\mathcal{S}} \overrightarrow{\boldsymbol{\mathcal{U}}}\_{\mathcal{S}} \boldsymbol{\varepsilon}\_{\mathcal{S}} \Big) = \nabla \cdot \Big( a\_{\mathcal{S}} \frac{\mu\_{t,\mathcal{S}}}{\sigma\_{\mathcal{E}}} \nabla \boldsymbol{\varepsilon}\_{\mathcal{S}} \Big) + a\_{\mathcal{S}} \frac{\mathcal{E}\_{\mathcal{S}}}{k\_{\mathcal{S}}} \Big( \mathbf{C}\_{1\varepsilon} \mathbf{G}\_{k,\mathcal{S}} - \mathbf{C}\_{2\varepsilon} \rho\_{\mathcal{S}} \boldsymbol{\varepsilon}\_{\mathcal{S}} \Big) \\ + a\_{\mathcal{S}} \rho\_{\mathcal{S}} \Pi\_{\varepsilon\_{\mathcal{S}}} \Big) \end{split} \tag{30}$$

In the above equations, *<sup>g</sup> <sup>k</sup>* and *<sup>g</sup>* represent the influence of the dispersed phase on the continuous phase and take the following forms:

$$\Pi\_{k\_{\mathcal{g}}} = \sum\_{p=1}^{m} \frac{K\_{\mathcal{g}^{\rm s}}}{\alpha\_{\mathcal{g}} \rho\_{\mathcal{g}}} \left( k\_{\mathcal{g}^{\rm s}} - \mathcal{D}k\_{\mathcal{g}} + \overrightarrow{\mathcal{U}}\_{\mathcal{g}^{\rm s}} \cdot \overrightarrow{\mathcal{U}}\_{\mathcal{d}^{\rm s}} \right) \tag{31}$$

$$
\Pi\_{\mathcal{E}\_{\mathcal{E}}} = \mathsf{C}\_{3\mathcal{E}} \frac{\mathcal{E}\_{\mathcal{g}}}{k\_{\mathcal{g}}} \Pi\_{k\_{\mathcal{g}}} \tag{32}
$$

The drift velocity *Udr* is defined in Equation (33). This velocity results from turbulent fluctuations in the volume fraction. When multiplied by the interchange coefficient *Kgs* , it serves as a correction to the momentum exchange term for turbulent flows:

$$\overrightarrow{\mathbf{U}}\_{dr} = -\left(\frac{D\_s}{\sigma\_{\mathcal{S}^\circ} a\_s} \nabla \alpha\_s - \frac{D\_\mathcal{g}}{\sigma\_{\mathcal{S}^\circ} a\_\mathcal{g}} \nabla \alpha\_\mathcal{g}\right) \tag{33}$$

such that *D DD <sup>g</sup> s ts*, *<sup>g</sup>* for Tchen Theory of multiphase flow (FLUENT 6.3 User's guide)[23]. The generation of turbulence kinetic energy due to the mean velocity gradients *Gk*,*<sup>g</sup>* is computed from:

$$\mathbf{G}\_{k,\mathcal{S}} = \mu\_{t,\mathcal{S}} \left( \nabla \cdot \overrightarrow{\mathbf{U}}\_{\mathcal{S}} + \left( \nabla \cdot \overrightarrow{\mathbf{U}}\_{\mathcal{S}} \right)^{T} \right) \colon \nabla \overrightarrow{\mathbf{U}}\_{\mathcal{S}} \tag{34}$$

The turbulent viscosity *<sup>t</sup>*,*<sup>g</sup>* given in the above equation is written in terms of the turbulent kinetic energy of the gas phase as:

$$
\mu\_{t,\mathcal{S}} = \rho\_{\mathcal{S}} \mathbb{C}\_{\mu} \frac{k\_{\mathcal{S}}^2}{\sigma\_{\mathcal{S}}} \tag{35}
$$

The Reynolds stress tensor defined in Equation (13) for the continuous phase is based on the Boussinesq hypothesis[24] given by:

$$\overset{\blacksquare}{\tau}\_{\mathcal{S}} = -\frac{2}{3} \Big( \rho\_{\mathcal{S}} \boldsymbol{a}\_{\mathcal{S}} k\_{\mathcal{S}} + \rho\_{\mathcal{S}} \boldsymbol{a}\_{\mathcal{S}} \mu\_{t,\mathcal{S}} \nabla \cdot \overrightarrow{\boldsymbol{U}}\_{\mathcal{S}} \Big) \overset{\blacksquare}{\mathbf{I}} + \rho\_{\mathcal{S}} \boldsymbol{a}\_{\mathcal{S}} \mu\_{t,\mathcal{S}} \Big( \nabla \cdot \overrightarrow{\boldsymbol{U}}\_{\mathcal{S}} + \left( \nabla \cdot \overrightarrow{\boldsymbol{U}}\_{\mathcal{S}} \right)^{T} \Big) \tag{36}$$

#### **2.5.2 Dispersed phase turbulence equations**

Time and length scales that characterize the motion of solids are used to evaluate the dispersion coefficients, the correlation functions, and the turbulent kinetic energy of the particulate phase. The characteristic particle relaxation time connected with inertial effects acting on a particulate phase is defined as:

$$\tau\_{F,sg} = \rho\_g \alpha\_s K\_{gs}^{-1} \left(\frac{\rho\_s}{\rho\_g} + \mathbb{C}\_V\right) \tag{37}$$

The Lagrangian integral timescale calculated along particle trajectories is defined as:

$$\tau\_{t, \text{sg}} = \frac{\tau\_{t, \text{g}}}{\sqrt{\left(1 + C\_{\beta} \xi^{2}\right)}}\tag{38}$$

where

94 Mass Transfer in Chemical Engineering Processes

( )

 

*H O l ds <sup>s</sup> <sup>X</sup> ds H O l H O l* 

<sup>2</sup>

To describe the effects of turbulent fluctuations of velocities and scalar quantities in each

In the context of gas-solid models, three approaches can be applied (FLUENT 6.3 User's guide)[23]: (1) modeling turbulent quantities with the assumption that both phases form a mixture of density ratio close to unity (mixture turbulence model); (2) modeling the effect of the dispersed phase turbulence on the gas phase and vice versa (dispersed turbulence model); or (3) modeling the turbulent quantities in each phase independent of each other (turbulence model for each phase). In many industrial applications, the density of the solid particles is usually larger than that of the fluid surrounding it. Furthermore, modeling the turbulent quantities in each phase is not only complex, but also computationally expensive when large number of particles is present. A more desirable option would then be to model the turbulent effect of each phase on the other by incorporating *source* terms into the conservation equations. This model is highly applicable when there is one primary phase (the gas phase) and the rest are dispersed dilute secondary phases such that the influence of the primary phase turbulence is the dominant factor in the random motion of the secondary

*t g g*

 

 

*U C G C*

*k kUU*

 

*t g g gg gg g gg k g g gg g g k k*

,

*gs dr k gs g*

3 *<sup>g</sup> <sup>g</sup> g k g*

 

*C k*

 

*g g gg g g g g g g k g gg*

for the effect of dispersed phase turbulence on the continuous phase as shown below:

 

*k Uk k G*

*t k*

1 2 *<sup>g</sup> <sup>m</sup> gs*

*p g g K*

 

 

turbulence models should be used for cases with swirl and vortex shedding (*RANS*, *k*

**2.5 Turbulence model equations** 

**2.5.1 Continuous phase turbulence equations**  In the case of multiphase flows, the standard *k*

,

 

> 

continuous phase and take the following forms:

phase, the *k*

phase.

and

*t*

 

In the above equations, *<sup>g</sup> <sup>k</sup>* and *<sup>g</sup>*

2 2

() ()

(28)

model equations are modified to account

, *<sup>g</sup>*

1, 2

   

(30)

 

 

*g*

represent the influence of the dispersed phase on the

(32)

*g*

 

*g g*

(31)

(29)

).

 

multiphase turbulent model can be used for simpler geometries. Advanced

$$\xi = \frac{\left| \overrightarrow{\boldsymbol{L}} \right|\_{\text{sg}} \left| \boldsymbol{\tau}\_{\boldsymbol{t}, \text{g}} \right|}{L\_{\text{t}, \text{g}}} \tag{39}$$

Numerical Simulation of Pneumatic and Cyclonic Dryers Using Computational Fluid Dynamics 97

desirable to simplify the computational domain to reduce computational time and effort and to prevent divergence problems. For instance, if the model shows some symmetry as in the case of a circular geometry, it can be modeled along the plane of symmetry. However, for a possible CFD solution to exist, the computational domain has to be discretized into cells or elements with nodal points marking the boundaries of each cell and combining the physical

It is a common practice to check and test the quality of the mesh in the model simply because it has a pronounced influence on the accuracy of the numerical simulation and the time taken by a model to achieve convergence. Ultimately, seeking an optimum mesh that enhances the convergence criteria and reduces time and computational effort is recommended. A widely used criterion for an acceptable meshing technique is to maintain the ratio of each of the cell-side length within a set number (x/y, y/z, x/z < 3). In practice, and for most computational applications, local residual errors between consecutive iterations for the dependent variables are investigated. In the case of high residual values, it is then recommended to modify the model input or refine the mesh properties to minimize

The choice of meshing technique for a specific problem relies heavily on the geometry of the domain. Most CFD commercial packages utilize a compatible pre-processor for geometry creation and grid generation. For instance, FLUENT utilizes Gambit pre-processor. Two types of technique can be used in Gambit, a uniform distribution of the grid elements, or what can be referred to as structured grid; and a nonuniform distribution, or unstructured grid. For simple geometries that do not involve rounded edges, the trend would be to use structured grid as it would be easier to generate and faster to converge. It should be noted that the number of elements used for grid generation also plays a substantial role in simulation time and solution convergence. The finer the mesh, the longer the computational time, and the tendency for the solution to diverge become higher; nevertheless, the higher

Based on the above, one tends to believe that it might be wise to increase the number of elements indefinitely for better accuracy in the numerical predictions on the expense of computational effort. In practice, this is not always needed. The modeller should always bear in mind that an optimum mesh can be attained beyond which, changes in the

 In the following, two case studies are discussed. In each case, the computational domain is discretized differently according to what seemed to be an adequate mesh for the geometry

Let us consider a 4-m high vertical pipe for the pneumatic drying of sand particles and another 25-m high vertical pipe for the pneumatic drying of PVC particles. For both cases, the experimental data, physical and material properties were taken from Paixao and Rocha (1998)[25] for sand, and Baeyens et al. (1995)[26] for PVC as shown in Table 1. Both models were meshed and simulated in a three-dimensional configuration as shown in

In Figure 1, hot gas enters the computational domain vertically upward, fluidizes and dries the particles as they move along the length of the dryer. As the gas meets the particles, particles temperature increases until it reaches the wet bulb temperature at which surface

domain into one computational entity.

the solution accuracy.

under consideration.

Figures 1 and 2.

**Case 1** 

numerical predictions are negligible.

these errors in order to attain a converged solution.

and

$$C\_{\beta} = 1.8 - 1.35 \left(\cos \theta\right)^{2} \tag{40}$$

In Equation (40), is the angle between the mean particle velocity and the mean relative velocity. The constant term *CV* = 0.5 is an added mass coefficient (FLUENT 6.3 User's guide)[23].

The length scale of the turbulent eddies defined in Equation (39) is given by:

$$L\_{t,g} = \sqrt{\frac{3}{2}} \mathcal{C}\_{\mu} \frac{k\_g^{3/2}}{\varepsilon\_g} \tag{41}$$

The turbulence quantities for the particulate phase include

$$k\_s = k\_\mathcal{g} \left(\frac{b^2 + \eta\_{\mathcal{sg}}}{1 + \eta\_{\mathcal{sg}}}\right) \tag{42}$$

$$k\_{\rm sg} = \mathcal{Z}k\_{\rm g} \left(\frac{b + \eta\_{\rm sg}}{1 + \eta\_{\rm sg}}\right) \tag{43}$$

$$D\_{t,sg} = \frac{1}{3} k\_{sg} \tau\_{t,sg} \tag{44}$$

such that

$$b = \left(1 + \mathcal{C}\_V\right) \left(\frac{\rho\_s}{\rho\_\mathcal{g}} + \mathcal{C}\_V\right)^{-1} \tag{45}$$

$$
\eta\_{\rm sg} = \frac{\tau\_{t,\rm sg}}{\tau\_{\rm F,sg}} \tag{46}
$$

### **3. Grid generation**

The development of a CFD model involves several tasks that are equally important for a feasible solution to exist with certain accuracy and correctness. A reliable model can only be possible when correct boundary and initial conditions are implemented along with a meaningful description of the physical problem. Thus, the development of a CFD model should involve an accurate definition of the variables to be determined; choice of the mathematical equations and numerical methods, boundary and initial conditions; and applicable empirical correlations. In order to simulate the physical processes occurring in any well defined computational domain, governing and complimentary equations are solved numerically in an iterative scheme to resolve the coupling between the field variables. With the appropriate set of equations, the system can be described in two- and three-dimensional forms conforming to the actual shape of the system. In many cases, it is desirable to simplify the computational domain to reduce computational time and effort and to prevent divergence problems. For instance, if the model shows some symmetry as in the case of a circular geometry, it can be modeled along the plane of symmetry. However, for a possible CFD solution to exist, the computational domain has to be discretized into cells or elements with nodal points marking the boundaries of each cell and combining the physical domain into one computational entity.

It is a common practice to check and test the quality of the mesh in the model simply because it has a pronounced influence on the accuracy of the numerical simulation and the time taken by a model to achieve convergence. Ultimately, seeking an optimum mesh that enhances the convergence criteria and reduces time and computational effort is recommended. A widely used criterion for an acceptable meshing technique is to maintain the ratio of each of the cell-side length within a set number (x/y, y/z, x/z < 3). In practice, and for most computational applications, local residual errors between consecutive iterations for the dependent variables are investigated. In the case of high residual values, it is then recommended to modify the model input or refine the mesh properties to minimize these errors in order to attain a converged solution.

The choice of meshing technique for a specific problem relies heavily on the geometry of the domain. Most CFD commercial packages utilize a compatible pre-processor for geometry creation and grid generation. For instance, FLUENT utilizes Gambit pre-processor. Two types of technique can be used in Gambit, a uniform distribution of the grid elements, or what can be referred to as structured grid; and a nonuniform distribution, or unstructured grid. For simple geometries that do not involve rounded edges, the trend would be to use structured grid as it would be easier to generate and faster to converge. It should be noted that the number of elements used for grid generation also plays a substantial role in simulation time and solution convergence. The finer the mesh, the longer the computational time, and the tendency for the solution to diverge become higher; nevertheless, the higher the solution accuracy.

Based on the above, one tends to believe that it might be wise to increase the number of elements indefinitely for better accuracy in the numerical predictions on the expense of computational effort. In practice, this is not always needed. The modeller should always bear in mind that an optimum mesh can be attained beyond which, changes in the numerical predictions are negligible.

 In the following, two case studies are discussed. In each case, the computational domain is discretized differently according to what seemed to be an adequate mesh for the geometry under consideration.

### **Case 1**

96 Mass Transfer in Chemical Engineering Processes

<sup>2</sup> *C* 1.8 1.35 cos

velocity. The constant term *CV* = 0.5 is an added mass coefficient (FLUENT 6.3 User's

3 2 3/2

*g*

*sg*

*sg*

*sg*

 

*sg*

1

*g*

*k*

 

2 1

 

1

, , 1 3 *D k t sg sg t sg* 

1 *<sup>s</sup> V V g*

, , *t sg*

The development of a CFD model involves several tasks that are equally important for a feasible solution to exist with certain accuracy and correctness. A reliable model can only be possible when correct boundary and initial conditions are implemented along with a meaningful description of the physical problem. Thus, the development of a CFD model should involve an accurate definition of the variables to be determined; choice of the mathematical equations and numerical methods, boundary and initial conditions; and applicable empirical correlations. In order to simulate the physical processes occurring in any well defined computational domain, governing and complimentary equations are solved numerically in an iterative scheme to resolve the coupling between the field variables. With the appropriate set of equations, the system can be described in two- and three-dimensional forms conforming to the actual shape of the system. In many cases, it is

*F sg*

*b*

 

*b*

is the angle between the mean particle velocity and the mean relative

(40)

(42)

(43)

(44)

(45)

(46)

(41)

The length scale of the turbulent eddies defined in Equation (39) is given by:

The turbulence quantities for the particulate phase include

,

*t g*

*L C*

*s g*

2

*sg g*

*bC C* 

*sg*

*k k*

*k k*

and

In Equation (40),

guide)[23].

such that

**3. Grid generation** 

Let us consider a 4-m high vertical pipe for the pneumatic drying of sand particles and another 25-m high vertical pipe for the pneumatic drying of PVC particles. For both cases, the experimental data, physical and material properties were taken from Paixao and Rocha (1998)[25] for sand, and Baeyens et al. (1995)[26] for PVC as shown in Table 1. Both models were meshed and simulated in a three-dimensional configuration as shown in Figures 1 and 2.

In Figure 1, hot gas enters the computational domain vertically upward, fluidizes and dries the particles as they move along the length of the dryer. As the gas meets the particles, particles temperature increases until it reaches the wet bulb temperature at which surface

Numerical Simulation of Pneumatic and Cyclonic Dryers Using Computational Fluid Dynamics 99

evaporation starts to occur. At this stage, convective mass transfer dominates the drying of surface moisture of particles during their residence time in the dryer. Since pneumatic drying is characterized by short residence times on the order of 1-10 seconds, mostly convective heat- and mass transfer occur. However, since experimental data for pore moisture evaporation were also provided in the independent literature, moisture diffusion

The computational domain was discretized into hexahedral elements with unstructured mesh in the *x* and *z*-directions and nonuniform distribution in the *y*-direction. An optimized mesh with approximately 63 000 cells and 411 550 cells was applied for the sand and PVC models, respectively. The computational grid is shown in Figure 2. Grid generation was done in Gambit 4.6, a compatible pre-processor for FLUENT 6.3. A grid sensitivity study was performed on the large-scale riser using two types of grids, a coarse mesh with 160 800 elements, and finer mesh with 411 550 elements. All models were meshed based on hexahedral elements due to their superiority over other mesh types when oriented with the direction of the flow. Results obtained for the axial profiles of pressure and relative velocity yield a maximum of 15% difference between the predicted results up to 4.5 m above the dryer inlet; however, there was hardly any difference in the results at a greater length by changing the size of the grids. Therefore, the coarsest grid was used in all simulations.

In this case, let us consider a different geometry as shown in Figure 3. This model discusses the drying of sludge material and linked to an earlier work presented by Jamaleddine and

or the second stage of drying was also considered.

Fig. 2. Computational grid

**Case 2** 


Paixao and Rocha (1998)[25]

Table 1. Conditions used in the numerical model simulation

Fig. 1. (Left) Geometrical models; (middle) sand model; (right) PVC model

evaporation starts to occur. At this stage, convective mass transfer dominates the drying of surface moisture of particles during their residence time in the dryer. Since pneumatic drying is characterized by short residence times on the order of 1-10 seconds, mostly convective heat- and mass transfer occur. However, since experimental data for pore moisture evaporation were also provided in the independent literature, moisture diffusion or the second stage of drying was also considered.

### Fig. 2. Computational grid

The computational domain was discretized into hexahedral elements with unstructured mesh in the *x* and *z*-directions and nonuniform distribution in the *y*-direction. An optimized mesh with approximately 63 000 cells and 411 550 cells was applied for the sand and PVC models, respectively. The computational grid is shown in Figure 2. Grid generation was done in Gambit 4.6, a compatible pre-processor for FLUENT 6.3. A grid sensitivity study was performed on the large-scale riser using two types of grids, a coarse mesh with 160 800 elements, and finer mesh with 411 550 elements. All models were meshed based on hexahedral elements due to their superiority over other mesh types when oriented with the direction of the flow. Results obtained for the axial profiles of pressure and relative velocity yield a maximum of 15% difference between the predicted results up to 4.5 m above the dryer inlet; however, there was hardly any difference in the results at a greater length by changing the size of the grids. Therefore, the coarsest grid was used in all simulations.

### **Case 2**

98 Mass Transfer in Chemical Engineering Processes

**Particle Sand PVC**  Diameter (mm) 0.38 0.18 Density (g / cm3) 2.622 1.116 Specific Heat [J / (kg oC)] 799.70 980.0

Height (m) 4.0 25.0 Internal Diameter (cm) 5.25 125.0 Gas Flow rate, Wg (kg/s) 0.03947 10.52 Solids Flow rate, Ws (kg/s) 0.00474 1.51 Inlet Gas Temperature, Tg (oC) 109.4 126.0 Inlet Solids Temperature, Ts (oC) 39.9 - Inlet Gas Humidity, Yg (kg/kg) 0.0469 - Inlet Moisture Content of Particles, Xs (kg/kg) 0.0468 0.206

**Drying Tube** 

Paixao and Rocha (1998)[25]

Table 1. Conditions used in the numerical model simulation

Fig. 1. (Left) Geometrical models; (middle) sand model; (right) PVC model

In this case, let us consider a different geometry as shown in Figure 3. This model discusses the drying of sludge material and linked to an earlier work presented by Jamaleddine and

Numerical Simulation of Pneumatic and Cyclonic Dryers Using Computational Fluid Dynamics 101

The numerical analysis is based on a 3D, Eulerian multiphase CFD model provided by FLUENT/ANSYS R12.0. Physical and material properties for the sludge material are shown in Table 2. The computational domain was discretized into hexahedral elements with approximately 230 385 cells. This element type was chosen as it showed better accuracy between the numerical predictions and experimental data than tetrahedral elements as shown in Bunyawanichakul et al. (2006)[28]. The computational grid is shown in Figure 4.

The governing equations along with the complementary equations are solved using a pressure based solution algorithm provided by FLUENT 6.3. This algorithm solves for solution parameters using a segregated method in such a manner that the equations are solved sequentially and in a separate fashion. Briefly stated, the solution parameters are initially updated. The *x*-, *y*-, and *z*-components of velocity are then solved sequentially. The mass conservation is then enforced using the pressure correction equation (SIMPLE algorithm) to ensure consistency and convergence of solution equations. The governing equations are spatially discretized using second-order upwind scheme for greater accuracy and a first-order implicit for time. This allows for the calculation of quantities at cell faces using a Taylor series expansion of the cell-centered solution about the cell centroid. More details related to this can be found in Patankar[29], or FLUENT 6.3 User Guide (2006)[23]. SAND AND PVC MODELS: A modified *k*-*ε* turbulence model is used along with the standard wall function for both phases in the vicinity of the wall. To avoid solution divergence, small time steps on the order of 1 × 10-4 to 1 × 10-6 are adopted. Solution convergence is set to occur for cases where scaled residuals for all variables fall below 1 × 10-

3, except for the continuity equation (1 × 10-4) and the energy equation (1 × 10-6).

standard wall function for both phases in the vicinity of the wall. Bunyawanichakul et al.[28] validated their numerical predictions with experimental data by adopting tetrahedral mesh

turbulence model is used along with the

Grid generation was done in Gambit 4.6, a compatible pre-processor for FLUENT.

Fig. 4. Computational grid

**4. Numerical parameters – numerical solvers** 

SLUDGE MODEL: For this model, a RNG *k*

Ray (2010)[3] for the drying of sludge in a large-scale pneumatic dryer. Material properties for sludge are shown in Table 2. The geometrical model is a large-scale model of a design presented by Bunyawanichakul et al. (2006)[28]. The computational domain consists of an inlet pipe, three chambers in the cyclone, and an outlet. Two parallel baffles of conical shape with a hole or orifice at the bottom divide the dryer chambers. As the gas phase and the particulate phase (mixture) enter the cyclone dryer tangentially from the pneumatic dryer, they follow a swirling path as they travel from one chamber to another through the orifice opening. This configuration allows longer residence times for the sludge thus enhancing heat- and mass-transfer characteristics.


**\***Sludge properties are taken from Arlabosse et al. (2005) [27]

Table 2. Conditions used in the numerical model simulation

Fig. 3. Schematic of the pneumatic-cyclone dryer assembly

The numerical analysis is based on a 3D, Eulerian multiphase CFD model provided by FLUENT/ANSYS R12.0. Physical and material properties for the sludge material are shown in Table 2. The computational domain was discretized into hexahedral elements with approximately 230 385 cells. This element type was chosen as it showed better accuracy between the numerical predictions and experimental data than tetrahedral elements as shown in Bunyawanichakul et al. (2006)[28]. The computational grid is shown in Figure 4. Grid generation was done in Gambit 4.6, a compatible pre-processor for FLUENT.

Fig. 4. Computational grid

100 Mass Transfer in Chemical Engineering Processes

Ray (2010)[3] for the drying of sludge in a large-scale pneumatic dryer. Material properties for sludge are shown in Table 2. The geometrical model is a large-scale model of a design presented by Bunyawanichakul et al. (2006)[28]. The computational domain consists of an inlet pipe, three chambers in the cyclone, and an outlet. Two parallel baffles of conical shape with a hole or orifice at the bottom divide the dryer chambers. As the gas phase and the particulate phase (mixture) enter the cyclone dryer tangentially from the pneumatic dryer, they follow a swirling path as they travel from one chamber to another through the orifice opening. This configuration allows longer residence times for the sludge thus enhancing

**Particle Sludge \***  Diameter (mm) 0.18 Density (kg / m3) 998.0 Specific heat [J / (kg oC)] 4182.0 Thermal Conductivity [W / (m oC)] 0.6

 Height (m) 8.0 Internal diameter (m) 6.0

heat- and mass-transfer characteristics.

**\***Sludge properties are taken from Arlabosse et al. (2005) [27]

Table 2. Conditions used in the numerical model simulation

Fig. 3. Schematic of the pneumatic-cyclone dryer assembly

**Drying Tube** 
