1 Introduction

Fluidized beds are widely used in many industrial fields such as petroleum, chemical and energy due to their uniform mixing, large contact surfaces between phases and efficient heat transfer and mass transfer (Van Der Hoef et al. 2006). Many studies and applications have shown that circulating fluidized beds have become an important reactor for the combustion and utilization of petroleum coke (Chen and Lu 2007). The process of biomass conversion and utilization is also often carried out in a fluidized bed (Cui and Grace 2007). However, fluidization of biomass particles is a difficult task because the non-spherical biomass particles are generally large and have a low density. In the actual industrial processes, inert particles such as silica sand, alumina and calcite are often added to the biomass to assist its fluidization (Fotovat et al. 2013). Therefore, understanding mixing behaviors of binary mixtures in a fluidized bed has specific significance for the design and optimization of related industrial processes.

In recent years, many researchers have conducted experimental studies on the macroscopic flow and mixing behaviors of granular systems containing non-spherical particles. Zhang et al. (2009) experimentally studied the mixing and separation behaviors of biomass–sand binary mixture in a fluidized bed and analyzed the variation of flow pattern and mixing index. Fotovat et al. (2014) investigated effects of superficial gas velocity and the mixture composition on the mixing/segregation characteristics by analyzing the circulatory motion of the active tracer. Shao et al. (2016) investigated the mixing behaviors in binary-component and multi-component particulate systems and found that the density of non-spherical particles has a more significant influence on the mixing characteristics than particle size and shape. Boer et al. (2018) developed a digital image analysis technique to explore the particle distributions in the fluidized bed and obtained qualitative and quantitative experimental results which can be used for validation of numerical models concerning non-spherical particle mixing. In general, it is difficult to obtain microscopic flow characteristics through experimental methods, and numerical simulation can obtain more detailed information, which is a good complement to the experiment. Among the many simulation methods, the computational fluid dynamics–discrete element method (CFD–DEM) has attracted much attention due to the fact that the motion characteristics of each particle can be obtained. Tsuji et al. (1993) coupled DEM proposed by Cundall and Strack (1979) with CFD to describe the gas–solid two-phase flow in fluidized beds. Currently, CFD–DEM has been widely used to evaluate the flow characteristics of spherical particles in bubbling fluidized beds (Müller et al. 2008; Luo et al. 2014; Wang et al. 2015; Wang et al. 2016a, b; Liu and van Wachem 2019), spouted beds (Takeuchi et al. 2004; He et al. 2016) and spouted fluidized beds (Luo et al. 2015; Van Buijtenen et al. 2011; Tang et al. 2017; Wang et al. 2016c).

In previous studies, the assumption of spherical particles was common. However, in the actual industrial process, particles are mostly non-spherical particles, and particle shape has a great influence on the macroscopic flow characteristics (Hilton et al. 2010). Recently, some scholars have applied CFD–DEM to study non-spherical particle flow behaviors in a fluidized bed. Hilton et al. (2010) studied the flow characteristics of spherical particles, cubes and ellipsoidal particles in a fluidized bed and found that the non-spherical particles have larger drag force at the same superficial gas velocity, leading to earlier fluidization. Zhou et al. (2011) evaluated the flow characteristics of ellipsoidal particles and analyzed the effect of aspect ratio on particle orientation and force structure. Cai et al. (2012) investigated the orientation distribution of cylindrical particles in the riser and found that the orientation of cylindrical particles is affected by their position more than the shape and local gas velocity.

In the process of fluidization and utilization of biomass, mixtures comprising non-spherical particles are often involved. Therefore, it is necessary not only to study the flow characteristics of the monocomponent non-spherical particle system, but also to further investigate the mixing characteristics of the multi-component non-spherical particle system. Ren et al. (2013) investigated the mixing behavior of binary mixtures comprising corn-shaped particles and spheres in a spouted bed and pointed out that particles with closer sizes could be better mixed. Vollmari et al. (2015) evaluated the mixing behavior bidisperse mixtures of various shaped particles and found that the mixing was strongly influenced by the shape of the particles and the addition of plates can accelerate the mixing. Ma and Zhao (2018) investigated the fluidization of binary mixtures comprising spheres and rod-like particles and analyzed the effect of volume fraction of the rod-like particles on the minimum fluidization velocity and particle orientation.

Currently, there are few simulation studies on the mixing characteristics of binary mixtures containing non-spherical particles in a fluidized bed, and the understanding of microscopic characteristics such as granular temperature and particle orientation in the bed is not sufficient. Therefore, it is necessary to further investigate the mixing behavior and microscopic characteristics of binary mixtures containing non-spherical particles in a fluidized bed. In this study, CFD–DEM simulations with rolling resistance were carried out to evaluate the mixing behaviors of binary mixtures comprising spherocylindrical particles and spherical particles in a fluidized bed. The macroscale mixing behaviors and microscale hydrodynamic characteristics including granular temperature and particle orientation were discussed.

2 Model description

2.1 Gas phase

The gas phase is considered to be a continuous phase and modeled by the volume-averaged Navier–Stokes equations as shown below:

$$\frac{{\partial (\varepsilon_{\text{g}} \rho_{\text{g}} )}}{\partial t} + \nabla \cdot (\varepsilon_{\text{g}} \rho_{\text{g}} \varvec{u}_{\text{g}} ) = 0$$
(1)
$$\frac{{\partial (\varepsilon_{\text{g}} \rho_{\text{g}} \varvec{u}_{\text{g}} )}}{\partial t} + \nabla \cdot (\varepsilon_{\text{g}} \rho_{\text{g}} \varvec{u}_{\text{g}} \varvec{u}_{\text{g}} ) = - \varepsilon_{\text{g}} \nabla P - \varvec{S}_{\text{p}} - \nabla \cdot \left( {\varepsilon_{\text{g}}\varvec{\tau}_{\text{g}} } \right) + \varepsilon_{\text{g}} \rho_{\text{g}} {\mathbf{g}}$$
(2)

where εg is the porosity; t is the time, s; ug is the velocity of gas phase, m/s; ρg is the density of gas phase, kg/m3; Sp is the momentum exchange source term, kg/(m2 s2); τg is the viscous stress tensor, kg/(m s2).

$$\varvec{S}_{\text{p}} = - \frac{1}{{V_{\text{c}} }}\sum\limits_{i}^{{N_{\text{c}} }} {\varvec{F}_{{\text{d}i}} }$$
(3)
$$\varvec{\tau}_{\text{g}} = \frac{2}{3}\mu_{\text{g}} \left( {\nabla \cdot \varvec{u}_{\text{g}} } \right)\varvec{I} - \mu_{\text{g}} \left[ {\left( {\nabla \varvec{u}_{\text{g}} } \right) + \left( {\nabla \varvec{u}_{\text{g}} } \right)^{\text{T}} } \right]$$
(4)

where Vc is the characteristic fluid volume, m3; Nc is the number of particles in the characteristic fluid volume; Fd is the drag force, N; μg is the gas shear viscosity, Pa s; I is the unit tensor.

2.2 Solid phase

The solid phase is treated to be a discrete phase, and the motion of each individual particle is directly solved by Newton’s second law. The translational motion equation of a non-spherical particle is given as follows:

$$m_{\text{p}} \frac{{{\text{d}}^{2} \varvec{r}}}{{{\text{d}}t^{2} }} = - V_{\text{p}} \nabla P + m_{\text{p}} {\mathbf{g}} + \varvec{F}_{\text{c}} + \varvec{F}_{\text{d}}$$
(5)

where mp is the mass of particle, kg; r is the position of particle center, m; vp is the translational velocity of particle, m/s; Fc is the contact force, N.

In this study, each non-spherical particle is modeled by the multi-sphere method (Favier et al. 1999). In the multi-sphere method, the contact force and torque acting on a non-spherical particle are obtained by calculating the sum of contact forces and torques on its spherical elements. Figure 1 shows the schematic diagram of transferring the force acting on the elemental sphere to the centroid of a non-spherical particle. To reduce the computational cost, a two-step method (Wang et al. 2016b) was used for contact detection.

Fig. 1
figure 1

Schematic diagram of transferring the force acting on the elemental sphere to the centroid of a non-spherical particle. a Forces, respectively, acting on each element sphere; b forces and torques acting on the center of each element sphere; c forces and torques acting on the center of a non-spherical particle

Full size image

The contact force acting on spherical element was calculated using a linear spring-dashpot (LSD) model and consists of a normal component Fn and a tangential component Ft:

$$\varvec{F}_{\text{n}} = - k_{\text{n}}\varvec{\delta}_{\text{n}} - \eta_{\text{n}} \varvec{v}_{\text{n}}$$
(6)
$$\varvec{F}_{\text{t}} = \left\{ {\begin{array}{*{20}l} { - k_{\text{t}}\varvec{\delta}_{\text{t}} - \eta_{\text{t}} \varvec{v}_{\text{t}} } \hfill &\quad {\left| {\varvec{F}_{\text{t}} } \right| \le \mu_{\text{f}} \left| {\varvec{F}_{\text{n}} } \right|} \hfill \\ { - \mu_{\text{f}} \left| {\varvec{F}_{\text{n}} } \right|\frac{{\varvec{v}_{\text{t}} }}{{\left| {\varvec{v}_{\text{t}} } \right|}}} \hfill &\quad {\left| {\varvec{F}_{\text{t}} } \right| > \mu_{\text{f}} \left| {\varvec{F}_{\text{n}} } \right|} \hfill \\ \end{array} } \right.$$
(7)

where k is the spring stiffness (N/m); δ is the elastic deformation (m); η is the damping coefficient (N s /m); v is the relative velocity (m/s); μf is the sliding friction coefficient. The subscripts n and t indicate the variables in normal direction and tangential direction, respectively.

Since the non-spherical particles inertia moment changes constantly during the rotation, the rotational motion of the non-spherical particles is more complicated than that of the spherical particles. To describe the rotational motion of non-spherical particles more conveniently, two coordinate systems are generally introduced: (a) a global coordinate system and (b) a local coordinate system. The angular momentum equation is solved in a local coordinate system, and the rotational variables are transformed between the local coordinate system and global coordinate system using the quaternion method (Wachem et al. 2015). As shown in Fig. 2, the direction of the coordinate axis of the local coordinate system is the direction of the principal axis of the particle’s inertia.

$$\varvec{I}_{\text{L}} \frac{{{\text{d}}\boldsymbol{\omega}_{\text{L}} }}{{{\text{d}}t}} + \boldsymbol{\omega}_{\text{L}} \times (\varvec{I}_{\text{L}} \cdot\varvec{\omega}_{\text{L}} ) = \varvec{T}_{\text{L}} = \varvec{T}_{\text{c}} + \varvec{T}_{\text{r}}$$
(8)

where IL is the inertia tensor (kg m2); ωL is the rotational velocity of particle in a local coordinate system (rad/s); TL is the torque in a local coordinate system (N m); Tc is the torque generated by contact force (N m); Tr is the torque generated by rolling friction (N m).

Fig. 2
figure 2

Schematic diagram of the local coordinate system

Full size image

The torque generated by rolling friction is attributed to the elastic hysteresis loss or viscous dissipation and slows down the relative rotations of particles (Zhou et al. 1999). Goniva et al. (2012) found that rolling friction can affect the velocity of spherical particles in spouted fluidized beds and simulation results can be improved significantly when applying a rolling friction model. In this study, the correlation proposed by Brilliantov and Pöschel (1998) was used to calculate the torque on contacting spherical elements generated by rolling friction as follows:

$$\varvec{T}_{\text{ri}} = - \mu r_{\text{p}} \left| {\varvec{F}_{\text{n}} } \right|\hat{\varvec{\omega }}_{i}$$
(9)

where μ is the rolling friction coefficient; rp is the radius of spherical element (m); |Fn| is the magnitude of normal contact force (N); ω is the rotational velocity of particle (rad/s).

2.3 Interphase exchange

Based on the drag force applied on a spherical particle, Di Felice (1994) proposed a drag correlation in multi-particle systems by experimental fitting.

$$\varvec{F}_{\text{d}} = \frac{1}{2}\rho_{\text{g}} C_{\text{D}} A_{ \bot } \left| {\varvec{u}_{\text{g}} \varvec{ - v}_{\text{p}} } \right|(\varvec{u}_{\text{g}} {\mathbf{ - }}\varvec{v}_{\text{p}} )\varepsilon_{\text{g}}^{2 - \chi }$$
(10)
$$\chi = 3.7 - 0.65\exp \left[ { - \frac{{(1.5 - \log_{10} {\text{Re}})^{2} }}{2}} \right]$$
(11)

Based on the experimental data obtained from Leith (1987), Ganser (1993), Tran-Cong et al. (2004) and a comprehensive numerical study, Hölzer and Sommerfeld (2008) established a simple and accurate drag model for non-spherical particles:

$$C_{\text{D}} = \frac{8}{\text{Re}}\frac{1}{{\sqrt {\varPhi_{\parallel} } }} + \frac{16}{\text{Re}}\frac{1}{\sqrt \varPhi } + \frac{3}{{\sqrt {\text{Re}} }}\frac{1}{{\varPhi^{{\frac{3}{4}}} }} + 0.42 \times 10^{{0.4( - \log \varPhi )^{0.2} }} \frac{1}{{\varPhi_{ \bot } }}$$
(12)

where \(\varPhi\) is the sphericity, \(\varPhi_{ \bot }\) is the crosswise sphericity, and \(\varPhi_{ \parallel}\) is the lengthwise sphericity. Re is defined as follows,

$${\text{Re}} = \frac{{\varepsilon_{\text{g}} d_{\text{p}} \left| {\varvec{u}_{\text{g}} - \varvec{v}_{\text{p}} } \right|\rho_{\text{g}} }}{{\mu_{\text{g}} }}$$
(13)

In this study, Di Felice (1994) correlation and Hölzer and Sommerfeld (2008) drag model were used to describe interphase momentum exchange. This method has been used in many studies in non-spherical granular systems (Hilton et al. 2010; Gan et al. 2017; Ma et al. 2017).

2.4 Particle granular temperature

The particle granular temperature is the turbulent energy of particles, indicating the velocity fluctuations of particles (Tartan and Gidaspow 2004; Jung et al. 2005). The particle translational granular temperature is defined as the mean of second moment of particle translational velocity fluctuation:

$$\left\langle {C_{i} C_{j} } \right\rangle = \frac{\text{1}}{n}\sum\limits_{{k = \text{1}}}^{n} {\left[ {v_{{{\text{p}},ik}} \left( {x,t} \right) - c_{i} \left( {x,t} \right)} \right]} \left[ {v_{{{\text{p}},jk}} \left( {x,t} \right) - c_{j} \left( {x,t} \right)} \right]$$
(14)

where i, j represent x, y or z directions, n is the number of particles in a unit volume, and vp(x,t) is the instantaneous translational velocity of particles. The average translational velocity ci(x,t) can be expressed as follows:

$$c_{i} \left( {x,t} \right) = \frac{\text{1}}{n}\sum\limits_{{k = \text{1}}}^{n} {v_{{{\text{p}},ik}} \left( {x,t} \right)}$$
(15)

Similarly, the particle rotational granular temperature (θp,rot(x,t)) can be obtained as follows:

$$\left\langle {C_{{\text{rot},i}} C_{{\text{rot},j}} } \right\rangle = \frac{\text{1}}{n}\sum\limits_{{k = \text{1}}}^{n} {\left[ {\omega_{{{\text{p}},ik}} \left( {x,t} \right) - c_{{\text{rot},i}} \left( {x,t} \right)} \right]} \left[ {\omega_{{{\text{p}},jk}} \left( {x,t} \right) - c_{{\text{rot},j}} \left( {x,t} \right)} \right]$$
(16)

where ωp(x,t) is the instantaneous rotational velocity of particles. The average rotational velocity crot,i(x,t) can be expressed as follows:

$$c_{{\text{rot},i}} \left( {x,t} \right) = \frac{\text{1}}{n}\sum\limits_{{k = \text{1}}}^{n} {\omega_{{{\text{p}},ik}} \left( {x,t} \right)}.$$
(17)

2.5 Initial and boundary conditions

As shown in Fig. 3, the dimension of the fluidized bed in the simulation was 150 mm, 20 mm and 400 mm in width, depth and height, respectively. At the bottom of fluidized bed, uniform airflow was injected through the air distribution plate. A pressure outlet boundary condition was used at the top of the bed, and a no-slip boundary condition was applied at the walls for the gas phase. The spring stiffness and damping coefficient significantly affect the maximum overlap of particle collisions. Kruggel-Emden et al. (2010) suggested a maximum overlap of 1% of the particle diameter to avoid alteration of the simulations. The spring stiffness was chosen so that the maximum overlap condition is fulfilled. The parameters applied in the simulation are summarized in Table 1.

Fig. 3
figure 3

Schematic geometry of the fluidized bed and the spherocylindrical particle

Full size image
Table 1 Parameters applied in the simulation
Full size table

The particles in the simulation included spherical particles (3 mm in diameter) and spherocylindrical particles modeled by a multi-sphere method (7.5 mm in long diameter, 3 mm in short diameter, as shown in Fig. 3). The volume fraction of the spherocylinders (Vol) in the binary mixtures in the fluidized bed was 100%, 75%, 50% and 25%, respectively. Composition of the granular systems in the simulation is summarized in Table 2. As shown in Fig. 4, at the initial moment (t = 0 s), the spherocylinders marked in blue are stacked above the spherical particles marked in red. In this study, the CFD time step is 5 × 10−5 s and the DEM time step is 5 × 10−6 s.

Table 2 Composition of the granular systems in the simulation
Full size table
Fig. 4
figure 4

Initial particle distribution in the fluidized bed (spherocylindrical particles marked in blue; spherical particles marked in red). a Vol = 100%, b Vol = 75%, c Vol = 50%, d Vol = 25%

Full size image

3 Results and discussion

3.1 Mixing behaviors

Figure 5 shows the instantaneous particle distribution at t = 10 s at Ug = 1.6 m/s and Ug = 2.0 m/s. It can be seen that under two superficial gas velocities, the spherocylindrical particles and the spherical particles are uniformly mixed at t = 10 s. In addition, due to the presence of the spherocylindrical particles, the large bubbles are not always in the middle of the bed and may appear in the vicinity of the left and right walls. It can also be observed that the bed expansion height is higher with superficial gas velocity of 2.0 m/s.

Fig. 5
figure 5

Snapshots of particle distribution at t = 10 s at Ug = 1.6 m/s and Ug = 2.0 m/s, aUg = 1.6 m/s (Vol = 100%, Vol = 75%, Vol = 50%, Vol = 25%), bUg = 2.0 m/s (Vol = 100%, Vol = 75%, Vol = 50%, Vol = 25%)

Full size image

The averaged translational velocities of the spherocylindrical particles and spherical particles during 0–10 s are shown in Fig. 6. The average particle translational velocities are calculated as follows:

Fig. 6
figure 6

The average translational velocities of the spherocylindrical particles and spherical particles during 0–10 s at aUg = 1.6 m/s and bUg = 2.0 m/s

Full size image
$$\left| {v_{\text{p}} } \right| = \frac{\text{1}}{N}\sum\limits_{{i = \text{1}}}^{N} {\sqrt {v_{{\text{p}x}}^{\text{2}} + v_{{\text{p}y}}^{\text{2}} + v_{{\text{p}z}}^{\text{2}} } }$$
(18)

It can be found that the instantaneous average translational velocities of the spherocylindrical particles at the initial stage are greater than that of the spherical particles, but after a few seconds, the average translational velocities of the two particles become very close. The time required for particle mixing can be roughly estimated by the moment when the average velocities of the two particles become very close. Particles mix faster when the superficial gas velocity is 2.0 m/s, compared with the superficial gas velocity of 1.6 m/s.

Moreover, Fig. 7 shows time-averaged translational velocities of the spherocylindrical particles and spherical particles during 5–10 s. It can be clearly seen that as the superficial gas velocity increases, time-averaged translational particle velocities in the bed increase. In addition, as the volume fraction of the spherocylindrical particles increases, time-averaged translational particle velocities in the bed also increase. This might be explained that non-spherical particles have larger drag at the same superficial gas velocity.

Fig. 7
figure 7

Time-averaged translational velocities of the spherocylindrical particles and spherical particles during 5–10 s at aUg = 1.6 m/s and bUg = 2.0 m/s

Full size image

The averaged rotational velocities of the spherocylindrical particles and spherical particles during 0–10 s are shown in Fig. 8. The average particle rotational velocities are calculated by:

Fig. 8
figure 8

The average rotational velocities of the spherocylindrical particles and spherical particles. aUg = 1.6 m/s, bUg = 2.0 m/s

Full size image
$$\left| {\omega_{\text{p}} } \right| = \frac{\text{1}}{N}\sum\limits_{{i = \text{1}}}^{N} {\sqrt {\omega_{{\text{p}x}}^{\text{2}} + \omega_{{\text{p}y}}^{\text{2}} + \omega_{{\text{p}z}}^{\text{2}} } }$$
(19)

It can be observed that even if the particles are uniformly mixed, the average rotational velocities of the spherical particles are larger than that of the spherocylindrical particles. This may be due to the fact that the spherocylindrical particles are in contact with more particles and the rotational resistance is larger.

In addition, Fig. 9 shows time-averaged rotational velocities of the spherocylindrical particles and spherical particles during 5–10 s. The results indicate that time-averaged rotational particle velocities in the bed increase with the increase in superficial gas velocity. As the volume fraction of the spherocylindrical particles increases, time-averaged rotational particle velocities in the bed also increase slightly.

Fig. 9
figure 9

Time-averaged rotational velocities of the spherocylindrical particles and spherical particles during 5–10 s at aUg = 1.6 m/s and bUg = 2.0 m/s

Full size image

To further analyze the particle mixing process, the Lacey mixing index (Lacey 1954) was used to quantitatively describe the mixing of binary particles:

$$M_{\text{Lacey}} = \frac{{\sigma_{0}^{2} - \sigma^{2} }}{{\sigma_{0}^{2} - \sigma_{r}^{2} }}$$
(20)

The variance defined as the square of the standard deviation can be calculated by the following formula:

$$\sigma^{2} = \frac{1}{N - 1}\sum\limits_{i = 1}^{N} {\left( {c_{i} - \bar{c}} \right)^{2} }$$
(21)

where N is the number of subcells; ci is the local concentration of sample particles; \(\bar{c}\) is the average concentration of sample particles; and \(\sigma_{0}^{2}\) and \(\sigma_{r}^{2}\) are the variances of the fully separated state and the completely mixed state, respectively.

The Lacey mixing index is 0 for the fully separated state, and the Lacey mixing index is 1 for the completely mixed state. The spherocylindrical particles are equally divided into equal volumes of spherical elements, and then, the Lacey mixing index is calculated using spherical elements and spherical particles. Figure 10 shows the Lacey mixing index during 0–10 s at the superficial gas velocity of 1.6 m/s and 2.0 m/s. It can be seen that the mixing index of both superficial gas velocities can reach a value close to 1, which indicates that the spherocylindrical particles and the spherical particles can achieve good mixing. However, corresponding mixing index can reach the stability more quickly when the superficial gas velocity is 2.0 m/s, which demonstrates that increasing the superficial gas velocity can accelerate the mixing of the particles. Moreover, when the superficial gas velocity is the same, as the volume fraction of the spherocylindrical particles increases, the particle mixing is slightly accelerated.

Fig. 10
figure 10

Variation of Lacey mixing index with time at Ug = 1.6 m/s and Ug = 2.0 m/s

Full size image

3.2 Particle granular temperature distribution

Figure 11 presents the time-averaged translational granular temperatures of spherocylindrical particles and spherical particles in z direction. The translational granular temperature is the translational turbulent energy of particles, indicating the translational velocity fluctuations of particles. It can be observed that the translational granular temperature of the spherical particles is slightly larger than that of the spherocylindrical particles. Moreover, as the volume fraction of the spherocylindrical particles increases, the translational granular temperature of the particles gradually increases, which means that the fluctuation of the translational velocity of the particles becomes more severe.

Fig. 11
figure 11

Translational granular temperature distribution in z direction (Ug = 2.0 m/s). a Spherocylindrical particles (Vol = 75%, Vol = 50%, Vol = 25%), b spherical particles (Vol = 75%, Vol = 50%, Vol = 25%)

Full size image

Similarly, the rotational granular temperature is the rotational turbulent energy of particles, indicating the rotational velocity fluctuations of particles. In a pseudo-2D fluidized bed, y direction is the primary direction for rotation of particles. Figure 12 shows the time-averaged rotational granular temperatures of spherocylindrical particles and spherical particles in y direction. It can be clearly seen that the rotational granular temperature of the spherical particles is much larger than that of the spherocylindrical particles. This indicates that the fluctuation of the rotational velocity of the spherical particles is much more severe than that of the spherocylindrical particles. In addition, as the volume fraction of the spherocylindrical particles increases, the rotational granular temperature of the particles gradually increases.

Fig. 12
figure 12

Rotational granular temperature distribution in y direction (Ug = 2.0 m/s), a spherocylindrical particles (Vol = 75%, Vol = 50%, Vol = 25%), b spherical particles (Vol = 75%, Vol = 50%, Vol = 25%)

Full size image

3.3 Particle orientation distribution

The orientation of the spherocylindrical particles significantly affects the fluidization and mixing of particles. In this paper, the orientation angle was defined on the X–O–Z plane and varied from 0° to 90° for quantification, similar to the definition of Zhou et al. (2011). 0° corresponds to the particles in horizontal direction (the long axis of the particle is parallel to the x-axis), and 90° represents the particles in vertical direction (the long axis of the particle is parallel to the z-axis). A schematic view of the orientation angle of spherocylindrical particles is shown in Fig. 13a. The orientation angle is divided into six ranges with an interval of 15° (i.e., 0°–15°, 15°–30°, 30°–45°, 45°–60°, 60°–75°, 75°–90°), and particles are counted at every angle range to calculate the percentage distribution of particle orientation.

Fig. 13
figure 13

Particle orientation distribution at Ug = 1.6 m/s and Ug = 2.0 m/s (Vol = 100%), a schematic of the orientation angle of the spherocylindrical particles, bUg = 1.6 m/s, cUg = 2.0 m/s

Full size image

Figure 13 shows the percentage distribution of particle orientation for Vol = 100% during 0–10 s at the superficial gas velocity of 1.6 m/s and 2.0 m/s. It can be seen that the spherocylindrical particles tend to be horizontally oriented at the initial packing state. When entering the fluidization stage, the tendency of the horizontal orientation of the spherocylindrical particles is weakened, and the tendency of the vertical orientation is increased. Moreover, the results indicate that as the superficial gas velocity increases, the tendency of the horizontal orientation of the rod-shaped particles decreases, and the particles tend to be more vertically oriented.

Figure 14 shows the percentage distribution of the orientation of the spherocylindrical particles at the superficial gas velocity of 2.0 m/s. It can be found that as the volume fraction of the spherocylindrical particles decreases, the percentage of the vertical orientation increases, that is, the addition of the spherical particles makes the spherocylindrical particles preferably distributed in a vertical orientation. There are two possible reasons for the vertical orientation of the spherocylindrical particles in a fluidized bed. On the one hand, it is related to the action of the drag force. When the spherocylindrical particles are horizontally oriented, the cross-sectional area in the vertical flow direction is large, resulting in a large drag force, which causes the particles to move in a direction in which the flow resistance is reduced. On the other hand, the addition of spherical particles increases the voids in the fluidized bed. The increase in the movement space of the spherocylindrical particles helps to reduce the interlocking phenomenon and promote the fluidization of the spherocylindrical particles. In conclusion, the combined effect of the two aspects makes the spherocylindrical particles in the binary granular system preferably oriented vertically.

Fig. 14
figure 14

Particle orientation distribution (Ug = 2.0 m/s), a Vol = 100%, b Vol = 75%, c Vol = 50%, d Vol = 25%

Full size image

In order to investigate the orientation distribution of the binary mixture in the vicinity of the walls, the time-average orientation of the spherocylindrical particles in each 10-mm region near the left- and right-side walls (marked in blue) is counted as shown in Fig. 15. It can be seen that the particles in the vicinity of the side walls are more vertically oriented due to the wall effect. The superficial gas velocity increases, and the vertical orientation trend is still obvious, indicating that the wall effect dominates the orientation behavior of the spherocylindrical particles in the vicinity of the walls.

Fig. 15
figure 15

Particle orientation distribution in the vicinity of the left and right walls, aUg = 1.6 m/s, bUg = 2.0 m/s

Full size image

4 Conclusion

In this study, CFD–DEM approach with rolling resistance was applied to evaluate the mixing behaviors of binary mixtures comprising spherocylindrical particles and spherical particles in a fluidized bed. The effects of the volume fraction of the spherocylindrical particles on the macroscopic mixing characteristics and microscopic flow behaviors including granular temperature and orientation distribution were discussed.

Increasing the superficial gas velocity accelerates the mixing process of the particles. Furthermore, when the superficial gas velocity is the same, as the volume fraction of the spherocylindrical particles increases, the particle mixing is slightly accelerated. Moreover, when well mixed, the average translational velocities of the spherocylindrical and spherical particles become very close, but the average rotational velocities of the spherical particles are larger.

In the fluidized bed, the translational granular temperature of the spherical particles is slightly larger than that of the spherocylindrical particles and the rotational granular temperature of the spherical particles is much larger than that of the spherocylindrical particles. Moreover, as the volume fraction of the spherocylindrical particles increases, the translational and rotational granular temperatures of the particles gradually increase.

In the binary granular system, as the volume fraction of the spherocylindrical particles decreases, the percentage of the vertical orientation increases, that is, the addition of the spherical particles makes the spherocylindrical particles preferably distributed in a vertical orientation. In addition, the wall effect dominates the orientation behavior of the spherocylindrical particles in the vicinity of the walls.