CFD Simulation of Particles Mixing and Segregation in a Rotary Coating Disk: Influence of Drag Models and Restitution Coefficient ()
1. Introduction
A rotary drum is typically a cylinder rotating about its central axis and offers appropriate fluid dynamics behaviour for the multiphase fluids contained within it. It is a critical piece of equipment used in processing granular materials in several engineering applications, such as kilns, mixers, dryers, grinding mills, and coating process reactors. The common use of rotary drums is also due to their ability to deal with large particle size distribution and significant differences in physical properties. However, the efficiency of rotating drums is highly dependent on fluid dynamic behaviour [1] [2]. The challenges of materials processing in rotary drums have been the subject of several research works [1] [2] [3] that have yielded important information on the dynamics of granular beds, particle movements, and heat transport through granular media, thereby providing knowledge of how to improve product quality and increase production. Depending on the rotational speed, filling degree, physical characteristics of granular materials, and drum geometry, rotating drums can exhibit six distinct flow regimes: slipping, slumping, rolling, cascading, cataracting, and centrifuging as shown in Figure 1 [4]. Each has its own flow behaviour, which adds to the complexity of its study.
The rolling and slumping regimes are the most employed in industrial applications, particularly in the rotary coating process. The best particle mixing is often achieved in the rolling mode [4] [5]. Based on the characteristics of the particle movement, the bed motion in the rolling regime can be further divided into two regions known as the active and passive layers [6]. Figure 2 depicts a schematic diagram of a rotary disk coater that operates in rolling mode. In the rolling mode, most of the particle mixing takes place in the active region, whereas mixing in the passive region is insignificant [7]. Particle mixing within the active layer determines the surface renewal rate that affects the rates of heat, mass transfer, chemical reactions [5] and controls the particle coefficient of variation in the coating process [8].
Although numerous researchers have done experimental studies on particle flow behaviour within rotating drums in response to various operational and design variables, this approach has inherent technical, time, and financial constraints for conducting comprehensive studies [9] [10] [11] [12].
On the other hand, the availability of high-capacity computers at an affordable price makes mathematical modelling a viable alternative to costly and time- consuming experimental investigations. Many of these mathematical modelling studies on rotary drums [1] [7] [13] [14] [15] have focused on dominated parameters such as drum rotational speed, filling degree, feed rate, and drum inclination angle that influence the rolling flow and the active layer thickness [2]. Besides the rotational speed and filling ratio pertaining to drum sizes, additional numerical parameters thought to have a significant influence on fluid bed behaviour and active-layer thickness include drag force, particle collision characteristics, and solid phase wall boundary conditions. Drag force, which reflects the interphase momentum transfer between the two phases in any multiphase flow system, has been widely studied in various other engineering applications [16] [17] [18].
Figure 1. The six regimes of granular bed motion in a rotating disk: (a) slumping, (b) surging, (c) rolling, (d) cascading, (e) cataracting, and (f) centrifuging [4].
Figure 2. Schematic diagram of a rotary coating disk operating in the rolling mode [1].
A few research groups have investigated the hydrodynamics of granular flow in rotating drums by studying the drag force that is exerted in other multiphase flow systems [4] [7]. Santos et al. [4] analyzed the effect of drag force on particle velocity profiles in a rotating drum using four different drag models: Syamlal-O’Brien, Gidaspow, Huilin et al. and Gibilaro et al. Karunarathne et al. [7] studied particle segregation in the transverse plane of a rotating cylinder using different particle-particle drag models: the Schiller-Naumann, Morsi-Alexander, and Syamlal-O’Brien Symmetric. All three drag models were capable of simulating particle segregation. The Schiller-Naumann and Morsi-Alexander models showed good agreement with the experimental results, whereas the Syamlal-O’Brien- Symmetric model showed some deviations. Some researchers [17] [19] [20] have also investigated the effect of restitution coefficient on solid-gas flow characteristics in fluidized beds. Mostafazadeh et al. [19] applied a multifluid Eulerian model using kinetic theory to investigate the effect of the restitution coefficient on the mixing and segregation of different particles sizes. As lower the values of the restitution coefficient decrease, the greater the loss of particle kinetic energy. As a result, both the mean velocity of the solid phase and the bed height decreased. Tagliaferri et al. [17] studied the mixing and segregation tendency of a binary mixture of particles, focusing on the effect of restitution coefficient and integration methods on bed dynamics. From their results, only the optimal value of 0.99 of the restitution coefficients affected the bubble volume fraction and jetsam concentration profiles. Taghizadeh et al. [20] investigated the effects of various factors such as rotation speed, restitution coefficient, and particle size on the hydrodynamic and particle granular temperature. The authors found that increasing the restitution coefficient raised the granular temperature at different depths of the bed, which affected the hydrodynamic behaviour of the bed.
It is essential to describe mixing and flow behavior accurately when it comes to coating process operations in the vital pharmaceutical, food, and fertilizer coating industries [21]. In controlled-release fertilizer (CRF) applications, uniform coating uniformity will have a positive effect on delayed release time to the soil and water retention. According to our previous discussion, there has been little attention paid to the influence of drag force and restitution coefficient on bed behavior, the thickness of the active-passive region, and the velocities of gas and solid phases in coating rotary drums. Therefore, the main purpose of this CFD study is to investigate the effects of both drag force and restitution coefficient on the hydrodynamics of a urea coating rotating disk, particularly in the active-passive region.
2. CFD Model & Simulation
This work was carried out applying the Eulerian-Eulerian multiphase flow modelling approach. The dynamic behavior of the gas-solid multiphase flow system can be described through a set of equations derived from continuity and momentum equations and the kinetic theory of granular flows.
2.1. Governing Equations in the Euler-Euler Method
Continuity Equation
Conservation of mass in the flow is represented by continuity equations for the gas phase and the solid phases,
(1)
(2)
where α is the volume fraction, ρ is density, and ν is velocity. The g and s refer to the gas and solids phases, respectively [21].
Momentum Equation
Momentum for the Gas Phase:
The influence of viscous, pressure, and gravity forces on the dynamics of the gas and solid particles are captured by the momentum equations as follows [21]:
(3)
(4)
where
is the slip velocity between the phases, Bgs denotes the drag force coefficient relevant to the phases g and s, P stands for the pressure, and
is the deviator effective stress tensor of fluid phase g gravity [21]:
Momentum for the Solid Phase:
(5)
Here, Ps, µs, ζs, and I are the solids pressure [Pa], the solids viscosity [Pa·s], the solids bulk viscosity [Pa·s] and the unit tensor, respectively. Ps represents the solid pressure (normal forces) [Pa] created due to particle-particle collisions in a flow due to presence of several solid phases [22].
(6)
ess is the particle-particle restitution coefficient between phase s and n [-], ds is the particle diameter [m], dss is the mean diameter of the particles in phase n and s [m], g0,ss and are the radial distribution function [-] and the granular temperature θs [m2/s2], respectively [22].
The bulk viscosity of the solids, in Equation (7), given by ζs
(7)
The solids shear viscosity in Equation (8) is given as,
(8)
The radial distribution g0,ss for solid phases can be expressed by Equations. (9) and (10) [16]:
(9)
(10)
2.2. Description of the Drag Models
This study assessed three drag models for the gas-solid interphase exchange coefficient: the Syamlal-O’Brien, Gidaspow, and Schiller-Naumann drag models.
Syamlal-O.Brien Drag Model
The Syamlal O.Brien drag model, given in Equation (11), was derived for a single spherical particle in a fluid, and modified with a relative velocity correlation (vr). The relative velocity correlation is the terminal settling velocity of a particle in a system divided by the terminal settling velocity of a single sphere [23]. The main idea behind this model is the assumption that the Archimedes number in a single particle and a multi-particle system is the same. The Archimedes number relates the gravitational forces to the viscous forces [23].
(11)
In this model, αg is the volume fraction of fluid or in this case gas, αs is the particle volume fraction, ρg is the fluid or gas density, ds is the particle diameter and |us-ug| is the absolute relative interracial velocity of the particles compared to the fluid. Since this model is derived for a single spherical particle, the drag factor CD is also a single particle model from Dalla Valle [24]. This CD is modified with the relative velocity correlation as shown in Equation (12).
(12)
Richardson and Zaki [25] analytical model based on experimental data was used for the relative velocity correlation vr. This model is given by [26], as shown in Equation (13).
(13)
(14)
(15)
The Reynolds number used in this model is the particle Reynolds number as given in Equation (16) [26].
(16)
Gidaspow Drag Model
The Gidaspow drag model is a combination of the Wen and Yu drag model Equation (17) and the Ergun model Equation (19) [23]. The Wen and Yu drag model uses a correlation extracted from Richardson and Zaki’s experimental data [25]. This correlation is valid when the internal forces are negligible, which implies that the viscous forces dominate flow behavior. The Ergun equation is derived for a dense bed and relates the drag to the pressure drop through porous media.
(17)
The drag factor CD is the drag factor for a spherical particle given by Equation (18) [23].
(18)
The Ergun equation is shown in Equation (19).
(19)
The Ergun equation is a combination of the Kozeny Carman and the Burke Plummer equations [27]. The Kozeny Carman is the first part of the Ergun equation and describes the viscous flow with low Reynolds number. The Burke Plummer equation describes the kinetic flow with high Reynolds number and is the second part of the Ergun equation. The constant is a shape factor for the particles. The combination of the two drag models (Equations (17) and (19)) in the Gidaspow drag model is shown in Equation (20) [23]
(20)
Schiller-Naumann Drag Model
The Schiller-Naumann model is the default method in the ANSYS Fluent, and it is acceptable for general use in all fluid-fluid multiphase simulations. When using the homogeneous population balance models, the interfacial area will be calculated directly from the population balance variables [21]. This model is used to simulate the drag between fluid phases in the multiphase flow. The drag function f is given as Equation (21) [21]. The relative Reynolds number Re for the primary phase g and secondary phase s is given as in Equation (21).
(21)
The drag coefficient CD is given as:
(22)
2.3. Geometry and Mesh Generation
In the Space Claim, a circle with 0.6 m diameter was created to represent the transverse plan of the rotating pan coater. A computational mesh created with 20,500 elements as shown in Figure 3.
2.4. Initial and Boundary Conditions
The main boundary condition of the transverse plane in the rotary pan coater is the relative motion between the bed material and the rotating wall. A no-slip condition was assumed, with the relative velocities of the gas and particles at the wall set to zero and particles subject to wall friction and gravity. Initially, all dispersed particles were located at the bottom of the pan coater. The bed height was 0.2 m at t = 0 s, with two granular solids of particles (urea and sand) and a filling limit of up to 40%. The physical properties of the gas and solid phases are given in Table 1.
Figure 3. Generated free triangular mesh for the 2D geometry.
Table 1. Physical properties of material and model parameters [21].
2.5. Solution Strategy and Convergence Criteria
The governing equations of the model were solved by the finite volume method. The computational domain was partitioned into a finite number of control volumes. The second order upwind discretization method was used for the discretization. Discretization is the method that converts the partial differential equations of the model to algebraic equations for numerical solution. This set of algebraic equations was solved using CFD-based code (Fluent 19.1). The “SIMPLE” algorithm was used to do the pressure-velocity coupling. Two solid phases were used to describe the binary mixture of particles. Volume fraction was discretized according to the QUICK scheme. Applying unsteady state conditions and Eularian model for laminar flow, the solver ran for 60 seconds with 10−3 s time step and 10−3 residual value for the convergence. Table 2 summarizes the simulation parameters and their corresponding values.
Table 2. Summary of the simulation parameters.
3. Results & Discussion
3.1. Restitution Coefficient Comparison
To achieve the statistical steady-state (pseudo steady-state) for solution, the time evolution of particles volume fraction in different region along the bed depth is monitored at varied restitution coefficients of 0.7, 0.8, 0.9, and 0.95. The restitution coefficient (ess) characterizes the energy dissipated during particles collision and is defined as in Equation (23). The lower the value of the restitution coefficient, the more kinetic energy is lost up on collision [19].
(23)
As seen in Figures 4-7, as time passes, the urea particles travel downward to the bed mid-region while the sand particles move upward to the bed top-region. The volume fractions vary rapidly at the start of the simulation (5 to 20 s) for the cases with ess = 0.7 and 0.8, but slowly after 30 s. On the other hand, with ess = 0.9 and 0.95 the simulation did not converge to a stationary state until the time exceeded 65 s. As a result, for the scenario when ess = 0.95, the time-averaged variables are computed between 65 and 75 s. The results show that as restitution coefficient values increase, so does the simulation time required for conforming to the stationary state. However, when the restitution coefficient increases, there is less dissipation of particle kinetic energy due to a more substantial elastic particle-particle collision. This may explain why, as the restitution coefficient increases, a longer simulation time is required before the mixing pattern reaches a stationary state. Furthermore, higher restitution coefficient value indicates the existence of a strong segregation phenomenon in the rotating process of the binary particle system.
Figures 4-7 also show that when the restitution coefficient increases, more urea particles accumulate at the mid-section of the bed indicating a stronger tendency for segregation. However, when the restitution coefficient increased to 0.95, the volume fraction value increased dramatically. A similar finding has been reported by Taghizadeh et al. [20] and Zhao et al. [31]. The restitution coefficients of 0.7 and 0.8 may be assigned to one class (A) and the 0.9 and 0.95 can be assigned to another class (B). The simulation results were also produced in the form of velocity vectors profiles of the granular phases, as shown in Figure 8. In Figure 8, irregular velocity vectors can be seen in the active region of the granular bed for the lower restitution coefficient values (class A), but the higher restitution coefficient values (class B) resulted in excellent rolling regime behaviour obtained.
3.2. Drag Models Comparison
The classical drag models available in Fluent 19.1 of Gidaspow, Syamlal-O’Brien, and Schiller-Naumann drag models were studied to find an optimal drag model for improved rotating bed hydrodynamics. To allow for complete segregation, the simulations were run for 60 s of real flow time. Figures 9-11 show the results obtained from the simulated three drag models. Considering the particle segregation, the simulations indicate that the Gidaspow and Syamlal-O’Brien drag
Figure 4. Urea volume fraction contour at ess = 0.7 at time step 5 s to 60 s.
Figure 5. Urea volume fraction contour at ess = 0.8 at time step 5 s to 60 s.
Figure 6. Urea volume fraction contour at ess = 0.90 at time step 5 s to 65 s.
Figure 7. Urea volume fraction contour at ess = 0.95 at time step 5 s to 75 s.
Figure 8. Sand velocity vectors at ess = 0.95, 90, 0.8 and 0.7.
Figure 9. Sand volume fraction and velocity vectors Gidaspow drag model.
Figure 10. Sand volume fraction and velocity vectors Syamlal-O’Brien drag model.
Figure 11. Sand Volume fraction and velocity vectors Schiller-Naumann drag model.
models can simulate the segregation contrary to Schiller-Naumann model. The results obtained with the Syamlal-O’Brien model shown in Figure 10 deviate from the results of the Gidaspow drag models. Although the model was able to simulate an acceptable rolling flow mode for the bed with larger active-region and less particle collision in the passive-region, it did not show the complete separation of the two particulate phases which was found by applying Gidaspow model (Figure 9). Even Though, Syamlal-O’Brien model shown is capable to predict the particle segregation in a rolling flow as in agree with [7] The minor deviation shown could be due to assumption that the Archimedes number is the same in a single particle and a multi-particle system for this model as predicted by [7] [23].
The Schiller-Naumann model (Figure 11), however, are not the proper models for this case, where there is no rolling flow mode was obtained, particles in the passive-region shows much collision when they move, and no completed segregation is expected to occur. This finding was in contrast with Karunarathne et al. [7] finding, were Schiller-Naumann model and the Morsi-Alexander are able to model the particle segregation phenomena showing a good agreement with the experimental results while the Syamlal-O’Brien-Symmetric model had some deviations. Were Santos et al. [4] demonstrated that, the influence of different drag models on particle velocity profile can be neglected in the case of a rotating drum operated in the rolling regime where there is no fluid entering or leaving the system. Disputes in the drag force results, reveal that further investigation is required.
3.3. Effect of Particle Size on Mixing
The effect of particle size on mixing of granular particles was studied using two particles of different size (100 μm and 1000 μm) and same densities. The density of the particles was the same as the density of urea solid granular. Figure12 shows the volume fraction of each phase of the mixture at 60 second. Variances in the particles size show a significant effect on each phase volume fraction. As predicated from Jha [32] study the particle size is a dominated parameter for segregation of small particles at different size ratio and mixing ratio. Where the smaller particles concentrated in the mid-section of the moving bed, while the larger particles sink towards the bottom of the rotating disk, similar results have been reported elsewhere [7] [32] [33].
3.4. Effect of Particles Density on Mixing
Mixing due to density variations was investigated using a constant particle size of 1000 μm. The study considered densities of 1300 kg/m3 and 1500 kg/m3. Figure 13 shows the volume fractions of the two particle types. From the figure, no
Figure 12. Volume fraction of particles with a diameter of (a) 1000 μm and (b) 100 μm in a mixture of particles with same density 1300 kg/m3 att = 60 s.
Figure 13. Volume fraction of particles with a density of (a) 1300 kg/m3 (b) 1500 kg/m3 in a mixture of particles with same diameter 1000 μm at t = 60 s.
major separation appears due to density differences whilst good mixing could be achieved under different particles density, this could be due to the insignificant different in the tested values. Moreover, a uniform mixing would be achieved at relative higher densities, similar to what was shown by Aissa et al. [12].
4. Conclusion
Commercial CFD software Fluent 19.1 was used to simulate a 2D granular rotating bed of rotary coating disk hydrodynamics characteristic, using an Eulerian-Eulerian CFD model. Drag models, restitution coefficient, particles size and density were investigated regarding segregation in a rolling flow. P-P restitution coefficient and drag model were shown to have a function in segregation and therefore in the hydrodynamic behaviour of the bed. Therefore, the overall aim of this study was to get the best available rotating bed hydrodynamics model and coefficient. The study discovered that more urea particles collect in the middle of the bed, which suggests the potential for segregation is much more obvious. The models show that Gidaspow and Syamlal-O’Brien drag models alone may replicate an appropriate rolling flow mode for the bed with bigger active regions and fewer particle collisions in the passive area. In the rolling phase, particle separation occurs owing to variations in particle size. Separation owing to condensation was less relevant when using various densities.
Acknowledgements
The research grants extended by Ministry of Higher Education (MOHE) Malaysia, and the research facilities provided by Universiti Teknologi PETRONAS, Malaysia, are highly acknowledged.