Abstract
Devices with impinging streams have been employed in various fields of chemical engineering, as a means of intensifying heat and mass transfer processes. The particle behavior in gas-particle two-phase impinging streams (GPISs), which is of essential importance for the research of transfer processes, was simulated by an Eulerian-Lagrangian approach in this paper. Collisional interaction of particles was taken into account by means of a modified direct simulation Monte Carlo (DSMC) method based on a Lagrangian approach and the modified Nanbu method. A quantitative agreement was obtained between the predicted results and the experimental data in the literature. The particle motion behavior and the distributions of particle concentration and particle collision positions were presented reasonably. The results indicate that the particle distribution in GPIS can be divided into three zones: particle-collision zone, particle-jetting zone, and particle-scattering zone. Particle collisions occur mainly in the particle-collision zone, which obviously results in a few particles penetrating into the opposite stream. The interparticle collision rate and the particle concentration reach their maximum values in the particle-collision zone, respectively. The maximum value of the particle concentration increases with the increasing inlet particle concentration according to a logarithmic function. The interparticle collision rate is directly proportional to the square of local particle concentration.
1. Introduction
The earliest emergence of impinging streams (ISs) can be traced back to the development and application of the Koppers-Totzek gasifier in 1953, and its scientific concept was proposed by Elperin [1] in the early 1960s. In an original gas-particle two-phase IS (GPIS), two gas-particle streams impinge against each other at high velocity. Particles from one stream penetrate into the other, which increases the relative velocity between the particles and the gas and prolongs the mean residence time of particles. This phenomenon yields the effect of intensifying the interphase heat and mass transfer in particulate systems [2]. Therefore, the GPIS devices have a broad application in many industrial processes, such as coal gasification, combustion, and drying [3]. The particle behavior in GPIS is of essential importance for research of the transfer processes; however, it is difficult to study by experimental means. This study will focus on the particle behavior in GPIS, such as particle motion and particle collision, from a viewpoint of mathematical modeling.
There are two classical models in simulation of gas-particle flow, the Eulerian-Eulerian two-fluid model which treats particles as a continuous phase and the Eulerian-Lagrangian particle trajectory model which tracks individual particles. For the former model, it is difficult to give different characteristics of particle motion, interparticle collisions, and interactions between gas and particles in different zones in GPIS, because this model is based on all kinds of hypotheses. Therefore, the particle trajectory model, which needs fewer assumptions, is used in this study to calculate the gas turbulence and the motion of individual particle. On the other hand, an important character in GPIS is the severe interparticle interaction which plays a significant role in the characteristics of the particle motion and concentration distribution. Therefore, the effect of collisions between particles should be taken into consideration in addition to the particle-fluid interaction in the numerical study of GPIS.
Many models have been reported for GPIS. In earlier studies, Kitron et al. [4] applied the particle Boltzmann transport equation [5, 6] including interparticle collisions to the study of GPIS. The Monte Carlo method was then used to solve this equation to obtain the velocity distribution of particle phase. However, the solution process of this method was too complex to be applied to the engineering practice. Guo et al. [7] and Ni et al. [8] used a Markov chain stochastic model to predict the residence time distribution of gas and particle in an opposed multiburner (OMB) gasifier. Other researchers [9–13] proposed various single-particle dynamics models to analyze the motion behavior of single particles in GPIS. These models neglected the interaction between particles and cannot give the particle collision effect.
The direct simulation Monte Carlo (DSMC) method, based on a Lagrangian approach developed in rarefied gas dynamics for handling collisions of a large number of gas molecules [14, 15], is an effective method of dealing with particle motion and collisions in dense gas-particle two-phase flow [16, 17]. The authors of this paper first applied the DSMC method, called as the traditional DSMC method in this paper, to GPIS to deal with the interparticle collision with acceptable computational cost [18]. Li et al. [19] then established a 3D model of the impinging zone of an OMB gasifier with DSMC method. This model reveals the concentration and the mean velocity profiles of particles in the impinging zone. Nevertheless, the particle flow has its own characteristics in GPIS, such as nonuniform distributions of particle concentration and particle collision positions, and a period of short time required for a particle to pass through the impinging zone. According to these characteristics of multiple spatial scales and multiple time scales, the DSMC method is modified to be more suitable to the numerical study of GPIS as in our earlier work [20]. This modified DSMC method is not restricted to flow field cells in the calculation of particle motion and collisions, and the particle time step in this method is adaptive and calculated using local particle parameters, such as the particle concentration and the velocities of particles nearby. In this work, this modified DSMC method is validated quantitatively using experimental results in the literature. Then, this method is applied to the numerical study of the particle behavior in GPIS. The distributions of particle concentration and particle collision positions are analyzed based on the calculation results. More specifically, the effect of the inlet particle mass flow rate will also be given. These results lay the foundations for further heat and mass transfer study in GPIS.
2. Mathematical Model
Real processes in practical engineering GPIS are very complex generally, which include multiphase flow, heat and mass transfer, and chemical reactions sometimes. To simplify the issue, the following assumptions are made in the present modeling: heat and mass transfer is ignored in flow process; no chemical reactions are considered; and particles are treated as rigid hard spheres with the same diameter.
2.1. Gas-Particle Two-Phase Flow
An Eulerian-Lagrangian approach is adopted to describe the gas-particle two-phase flow in GPIS as mentioned above. The interactions between gas phase and particle phase are determined by means of Newtonian third law. Therefore, four-way coupling between discrete phase and continuous phase is carried out. The key features of the mathematical model are briefly descried as follows.
2.1.1. Governing Equations for Gas Phase
The gas flow field in GPIS is generally a turbulent flow field with backflows, and it is calculated using the realizable model [21] in this simulation. In this model, the conservation equations for mass and momentum can be written as follows: where is the external body forces that arise from interaction with the particle phase, which results in the momentum exchanges between the gas and particle phases. And is calculated from the reverse action of particles. To close (2), the modeled transport equations for the turbulence kinetic energy, , and its dissipation rate, , are where ; ; and represents the generation of turbulence kinetic energy due to the mean velocity gradients. is the generation of turbulence kinetic energy due to buoyancy. represents the contribution of the fluctuating dilatation in compressible turbulence to the overall dissipation rate. The model constants are , , , and [22].
2.1.2. Particle Motion Equations
In GPIS, the spherical particles are entrained by gas flow. The translational motion of a particle in the gas phase is governed by Newton’s second law of motion and can be written as where is an additional force term used in the simulation, including virtual mass force, pressure gradient force, and Saffman’s lift force. In addition, represents the drag force of gas phase acting on the particle which is the most important force acting on the particle calculated using the following equation: The drag coefficient, , is determined by the following equation: where Re is the relative Reynolds number, which is defined as
2.2. Interparticle Collision
The modified DSMC method, based on a gridless approach and local time stepping with automatic adaption, is applied for modeling the interparticle collisions in GPIS here. The basic ideas of this method are as follows:(1)real particles are replaced by smaller number of sampled particles. And the trajectories of sampled particles are calculated. In this way, the computer memory and computation time can be saved significantly. This item is just the same as the traditional DSMC method;(2)sampled particle motion is decoupled into the movement and collision processes as the traditional DSMC method. Sampled particle movement obeys the single particle motion model, and the collision process follows the interparticle collision dynamics. In the modified DSMC method, the particle time step, that is, the time required for a particle in each movement, is adaptive, which is different from the traditional DSMC method;(3)collision pairs are found through the theory of collision probability instead of using trajectories as the traditional DSMC method. Nevertheless, the calculation of the collision probability between particles and the searching of collision pairs is based on the spherical region generated from the local particle parameters rather than the flow field cells.
To update the particle time step, the local collision mean-free path, that is, the average distance a particle travels between collisions with other moving particles, is calculated firstly using the following equation: where is the velocity magnitude of the particle being tracked. is the collision frequency, that is, the number of collisions occurred in unit time, which is obtained according to the local distribution of particles. To satisfy the “principle of uncoupling” in the modified DSMC method, the particle time step is given by where , , and are components of particle velocity . is the time step for gas phase.
In the modified DSMC method, the collision pairs are found through the spherical region, that is, searching scope, with the particle tracked being considered as center and as its radius given below where is the maximum value of relative velocity between the particle being tracked and the sampled particles nearby. Then, the probability of collision between particle being tracked and particle during a time step is calculated by where and are the diameters of particle and particle , respectively. is the relative velocity vector between particles and . is the real particle number represented by the sampled particle . is the volume of the searching scope and equal to .
The modified Nanbu method [23] is used to search a candidate collision partner , which probably collide with particle in the searching scope during a time step , and decide whether the collision occurs. First, a random number with a uniform distribution from zero to unity is extracted from a generator. The candidate collision partner is selected according to where is defined as the integer part of and is the total number of all sampled particles in the searching scope. If the relation is satisfied, particle would collide with particle during the time step , and the velocities of particles and are replaced with the postcollision velocities, but without changing their positions. The post-velocities of the two particles are as follows: when when where is the normal unit vector directed from particle to particle at the moment on contact. is the unit vector in the tangential direction. is the tangential component of the relative velocity.
2.3. Numerical Solution
The numerical solutions of the gas-particle two-phase flow, including the governing equations for gas continuous phase, motion equations for sampled particles, and the momentum exchanges between two phases, are obtained using the CFD code FLUENT. The calculations of particle time step and interparticle collision process are performed using self-written program in VC language, which is taken as self-defined function (UDF) embedded in the calculation of gas-particle two-phase flow in GPIS. To obtain the rapid convergence of the calculation of gas-solid two-phase flow, the gas flow field is calculated to be convergent first, and then the particles are introduced. When a particle reaches the outlet, the particle tracking is terminated. If the number of particles in GPIS does not change greatly, the calculation is assumed to be convergent.
3. Model Validation
To validate the established model, two experimental cases (Case 1 and Case 2) in the literature [24] are calculated first. Figure 1 shows the schematic diagram of the experimental apparatus with two opposed nozzles. The inner diameter of either nozzle was 8 mm. The distance between two nozzles was 140 mm. Particles were conveyed by air into the nozzles. The parameters of the gas and particle phases injecting from nozzles for two cases are also summarized in Table 1.
The coordinate system is shown in Figure 1. The origin of the coordinate system is placed on the midpoint of the line connecting the two nozzles. The axial direction -axis is on the horizontal line through nozzle axes and defined as positive rightward. The -axis is called the axial line of the GPIS device. Figure 2 shows the meshed geometry for the calculation domain, and the domain is meshed with about 80000 hexahedral cells. The mesh refinement in the zone between two nozzles is carried out to raise the accuracy of gas phase solutions because of the big variation gradients of various variables in this zone. The minimum mesh size is about 1.3 mm.
Figure 3 shows the comparison between experimental and simulated dimensionless particle concentration profiles along the axial line, in which is the particle concentration at the outlet of nozzle. It can be seen that the particle concentration decreases sharply with the increasing distance from the outlet of nozzle along the axial line because of dispersion, and it increases obviously near the center because of the particle collisions. The predicted results of the particle concentration are in reasonable agreement with the experimental data, which indicates the validity of the present method.
(a) Case 1
(b) Case 2
4. Results and Discussion
Figure 4 shows the sketch of a laboratory-scale GPIS device with two opposed nozzles [25]. The impinging chamber is made of stainless steel of 0.5 m in diameter. The distance between the two nozzles is set to 200 mm; both nozzles have a diameter of 42 mm. The outlet diameter of the impinging chamber is relatively small to collect particles. The following part of work is to model the particle behavior in the impinging zone. To reduce the computational cost, the calculation domain is located near the nozzles framed by the black dash-line rectangle as shown in Figure 4. The calculation domain is meshed with hexahedral cells. The mesh refinement in the zone near the height of nozzles is carried out to raise the accuracy of gas phase solutions. In order to obtain grid-independent solutions, several mesh sizes, that is, 62740, 142394, and 274836, were tested. The results of the particle concentration along the axial line computed using 142394 and 274836 meshes were closer to each other; the maximum differences in the predicted results were indeed less than 10%. On the other hand, finer meshes required excessive computational time; therefore, mesh size of 142394 was selected in this study (see Figure 4).
The parameters used in the simulation are summarized in Table 2 with reference to the parameters used practically. The boundary conditions are “velocity-inlet” and “outflow” for inlet and outlet, respectively. All walls are treated as nonslip boundaries in the lateral direction for the gas phase with standard wall function. To obtain the rapid convergence of the calculation of gas-solid two-phase flow, the gas flow field is calculated to be convergent first, and then the particles are introduced. When a particle reaches the outlet, the particle tracking is terminated. If the number of particles in GPIS does not change greatly, the calculation is assumed to be convergent.
The coordinate system is shown in Figure 4. The origin of the coordinate system is placed on the midpoint of the line connecting two nozzles. The rightward -axis on the nozzle axes is called the axial line of the GPIS device. In addition, coordinate surface is called impinging plane.
4.1. Particle Motion Behavior
Figure 5 shows the particle position distribution at different moments for the case with kg/s in which the colormap denotes the particle residence time, and Figure 6 shows the velocity distribution of particles between two nozzles. From those two figures, it is found that the particle velocity remains about the same before the time point 0.01 s at which two particle flows start to meet each other near the origin of the coordinate system. This is mainly because that the collision between particles is caused by difference in particle velocities. Before this time point, particles move towards the impinging plane along the axial direction. The velocity difference between neighboring particles is small, which results in the small interparticle collision probability. In addition, the change in particle velocities is not obvious due to the small impacting impulse if two particles collide with each other.
(a) 0.006 s
(b) 0.01 s
(c) 0.02 s
(a) Particle velocity component on the axial line
(b) Particle velocity component in the impinging plane
When two gas-particle flows with a great number of particles meet each other at 0.01 s in the area around the origin of the coordinate system, particles from the opposite streams collide with each other violently. This area is called impingement zone generally. The collisions in this zone increase the velocity components in the impinging plane and decrease the axial velocity components of particles. Therefore, particles cannot pass through this area and penetrate into the opposite stream easily. These particles spread out from this area after several collisions as the particle position distribution at 0.02 s shown in Figure 5(c).
4.2. Particle Concentration Distribution
Figure 7 shows the residence time distribution of particles between two nozzles. It is found that the violent particle collisions increase the particle residence time in the impingement zone. In other words, the violent particle collision makes the particles accumulate in this zone, which results in the highest particle concentration near the origin as shown in Figure 8. And the high particle concentration strengthens the particle collision in turn. After particles leave the impingement zone in all directions, the particle concentration decreases drastically because of the increasing volume occupied by the particles.
| (a) Horizontal plane through nozzles () |
| (b) Impinging plane () |
Figure 9 shows the particle concentration profiles along the axial line for different inlet particle mass flow rates. It can be seen that particle concentration along the axial line increases with the increasing inlet particle mass flow rate reasonably. For a certain case, the particles injected into the device move towards the impinging plane along the axial line and the collisions with neighboring particles obviously can not change the particle paths as mentioned above. Therefore, the particle concentration does not change significantly along the axial line in the areas between either the nozzle or the impingement zone. In the impingement zone, the particle concentration increases sharply and reaches a maximum value near the origin due to the particle collisions. The relation between the maximal particle concentration in the impingement zone and the inlet particle concentration is drawn in Figure 10. It can be seen that the maximum value of particle concentration increases sharply with the increasing inlet particle concentration for lower . The relative particle concentration, ratio of and , also increases sharply at first. With the further increase in , the increase in with tends to slow down, and the ratio of and decreases obviously. The change of with can be fitted by a logarithmic function with a correlation coefficient 0.9860 as follows: where , , and . The ratio of with reaches the maximum value as kg/m3.
4.3. Particle Collision Position Distribution
Figure 11 shows the distribution of particle collision positions in 0.01 s for the case with kg/s, where one red point represents one collision occurring there. In the 0.01 s, there are 184059 collisions between sampled particles in total. According to the distribution of particle collision positions, the particle distribution can be divided into three zones: particle-collision zone, particle-jetting zone, and particle-scattering zone. The statistical result shows that particle collisions occur mainly in the particle-collision zone in which the number of collisions accounts for about 90% of the total. The particle-collision zone is an ellipsoid whose center is located at the origin of coordinates. The shorter axis on -axis and the longer axis on -axis of the ellipsoid are about 0.075 and 0.1 for the case with kg/s, respectively. In this zone, particles from the opposite streams collide with each other violently, which causes two results. First, particle collisions make the particles accumulate in this zone, which results in the highest particle concentration shown in Figure 8. Second, particle collisions make the particles spread out from this zone to the whole device. Between the particle-collision zone and two nozzles, the two long and narrow zones are called particle-jetting zone. Although the particle concentration in this zone is high, the smaller collision probability between particles caused by small velocity difference results in the smaller number of collisions in these two zones. The broad zone outside the particle-collision zone and the particle-jetting zones is called the particle-scattering zone. Particles spreading out from the particle-collision zone scatter quickly in this zone. The collision probability between particles in this zone is small because of the low particle concentration shown in Figure 8. Therefore, there are only a few collisions occurring near the particle-collision zone and few collisions occur in the area away from the particle-collision zone.
Figure 12 shows the profiles of interparticle collision rate along the axial line for different inlet particle mass flow rates, where the interparticle collision rate is defined as the number of interparticle collisions per unit time per unit volume. It is found that the profiles of interparticle collision rate along the axial line for different inlet particle mass flow rates are similar to those of particle concentration shown in Figure 9, and the maximum value of the interparticle collision rate is reached in the particle-collision zone due to the violent collisions. This phenomenon indicates that the interparticle collision rate has direct relation with the particle concentration, and the violent particle collisions occur in the particle-collision zone with high particle concentration. Figure 13 shows the relation between the average interparticle collision rate at the nozzle outlet and the local particle concentration for different cases. The fitting curve is as follows: where is a constant and equal to . The correlation coefficient of this fitting curve is 0.9966. From this formula, it can be seen that the interparticle collision rate is directly proportional to the square of particle concentration, which is in accordance with the statement in the literature [26].
5. Conclusions
The modified direct simulation Monte Carlo (DSMC) method is applied to predict the particle behavior in GPIS in this work. The main results are as follows:(1)according to the distribution of particle collision positions, the particle distribution can be divided into three zones: particle-collision zone, particle-jetting zone, and particle-scattering zone. Particle collisions occur mainly in the particle-collision zone, and the distribution of particle collision positions is very similar to the particle concentration distribution;(2)for GPIS, a few particles penetrate into the opposite stream due to the violent particle collisions and almost all particles spread out from the particle-collision zone;(3)the particle concentration and the interparticle collision rate along the axial line have a similar variation trend. They increase sharply and reach their maximum values near the origin in the particle-collision zone. The maximum value of the particle concentration increases with the increasing inlet particle concentration according to a logarithmic function. The interparticle collision rate is directly proportional to the square of local particle concentration.
Nomenclature
| : | Constants in fitting curve |
| : | Particle concentration, kg/ |
| : | Constants in turbulence model |
| : | Particle concentration on the axial line, kg/ |
| : | Drag coefficient |
| : | Diameter of particle, m |
| : | Diameter of nozzle, mm |
| : | Coefficient of restitution |
| : | Force vector, N |
| : | Acceleration vector of gravity |
| : | Relative velocity vector, m/s |
| : | Turbulence kinetic energy, / |
| : | Collision mean-free path, m |
| : | Distance between nozzles, mm |
| : | Mass of particle, kg |
| : | The number of the sampled particles |
| : | Interparticle collision rate, 1/-s |
| : | Normal unit vector |
| : | Number of real particles represented by a sampled particle |
| : | Static pressure, Pa |
| : | Collision probability |
| : | Random number on interval |
| : | Relative Reynolds number |
| : | Radius of searching scope, m |
| : | Simulation time, s |
| : | Tangential unit vector |
| : | Time step, s |
| : | Gas velocity vector, m/s |
| : | Gas velocity magnitude, m/s |
| : | Volume of searching scope, |
| : | Particle velocity vector, m/s |
| : | Particle velocity magnitude, m/s |
| : | Inlet particle mass flow rate for single nozzle, kg/s. |
| : | Turbulent dissipation rate |
| : | Dynamic viscosity, kg/m-s |
| : | Turbulent viscosity, kg/m-s |
| : | Friction coefficient of Coulomb’s friction law |
| : | Gas density, kg/ |
| : | Inverse effective turbulent Prandtl numbers for and |
| : | Stress tensor |
| : | Collision frequency. |
| (0): | Before the collision. |
| 0: | At nozzle outlet |
| : | Tangential direction |
| : | Drag force |
| : | Gas phase |
| : | Sampled particle No. |
| : | Impinging plane |
| : | Turbulent kinetic energy |
| max: | Maximum value |
| : | Particle |
| : | Relative velocity |
| : | Vector direction |
| : | Turbulent dissipation rate. |
Acknowledgments
The financial supports from the National Natural Science Foundation of China (no. 50976024), the National Basic Research Program of China (no. 2010CB227002), the Natural Science Foundation of Jiangsu Province (no. BK20130520), and the Senior Expert Fund of Jiangsu University (no. 11JDG152) are greatly appreciated.