Generalized Boltzmann kinetic theory for EMMS-based two-fluid model

Chemical Engineering Science 156 (2016) 44–55

Contents lists available at ScienceDirect

Chemical Engineering Science journal homepage: www.elsevier.com/locate/ces

Generalized Boltzmann kinetic theory for EMMS-based two-fluid model Bidan Zhao a,b, Shuyue Li a,c, Junwu Wang a,n a

State Key Laboratory of Multiphase Complex Systems, Institute of Process Engineering, Chinese Academy of Sciences, Beijing 100190, PR China University of Chinese Academy of Sciences, Beijing 100049, PR China c State Key Laboratory of Heavy Oil Processing, China University of Petroleum, Beijing 102249, PR China b

H I G H L I G H T S








Kinetic theory is used to derive the governing equation of gas–solid flow. Molecular and pseudo-Reynolds stresses of gas phase can be obtained simultaneously. Boltzmann equation for the variation of volume, density and velocity of clusters is derived. Statistical mechanics foundation of EMMS-based two-fluid model is developed. CFD simulations using EMMS-based two-fluid model are performed.

art ic l e i nf o

a b s t r a c t

Article history: Received 29 March 2016 Received in revised form 5 September 2016 Accepted 8 September 2016 Available online 9 September 2016

It has long been recognized that the solid particles in circulating fluidized bed risers are distributed heterogeneously in the form of clusters. In response to this fundamental phenomenon, an EMMS-based two-fluid model has been developed recently from the viewpoint of continuum mechanics, however, its microscopic foundation remains unknown. In this study, the statistical mechanics foundation of EMMSbased two-fluid model was presented using generalized Boltzmann kinetic theory. With respect to the gas phase, a new method was developed by considering the fluctuations at different scales simultaneously, with which we can for the first time derive the correct governing equations of gas phase via kinetic theory, in the sense that both the molecular stress and the Reynolds (or pseudo-Reynolds) stress can be obtained simultaneously, whereas all previous kinetic theory analyses failed to predict the appearance of Reynolds (or pseudo-Reynolds) stress in the momentum conservation equation of gas phase due to the assumption of uniform structure, although it is physically always existent no matter how small the Reynolds number is. In case of particle phase, the generalized Boltzmann equation considering the spatio-temporal variation of the volume, density and velocity of clusters was firstly derived, a set of macroscopic transport equations was then derived in different phase spaces. It was shown that the governing equations of dense phase in the EMMS-based two-fluid model derived from continuum mechanics viewpoint corresponds to the macroscopic transport equations at (r, t ) space. Therefore, present study launches a solid microscopic foundation of EMMS-based two-fluid model. Finally, CFD simulations have been carried out to validate EMMS-based two-fluid model and to study the effect of gas phase pseudo-turbulence. & 2016 Elsevier Ltd. All rights reserved.

Keywords: Fluidization Multiphase flow Pseudo-Reynolds stress Mesoscale structures Cluster Kinetic theory

1. Introduction Circulating fluidized beds (CFBs) have been widely used in industry (Reh, 1971; Kunii and Levenspiel, 1991; Reh, 2003), such as coal combustion and fluid catalytic cracking (FCC). A variety of n

Corresponding author. E-mail address: [email protected] (J. Wang).

http://dx.doi.org/10.1016/j.ces.2016.09.012 0009-2509/& 2016 Elsevier Ltd. All rights reserved.

experimental studies on the hydrodynamics of gas–solid flow in CFB risers have been conducted (Reh, 1971; Yerushalmi et al., 1976; Berruti et al., 1995; Lim et al., 1995), it was undoubtedly concluded that the solid particles in CFB risers are distributed heterogeneously, which takes the form of clusters (Li et al., 1991; Zou et al., 1994; Soong et al., 1994; Li et al., 1995; Lackermeier et al., 2001; Cocco et al., 2010). Furthermore, the size of clusters, the solid concentration inside clusters (or cluster density) and the

B. Zhao et al. / Chemical Engineering Science 156 (2016) 44–55

velocity of clusters are varied spatio-temporally (Soong et al., 1994; Sharma et al., 2000; Harris et al., 2002; Manyele et al., 2002; Wang et al., 2008), and those particle clustering structures play a critical role in determining the hydrodynamics of gas–solid flow in CFB risers (Li et al., 1988; Li and Kwauk, 1994; Li et al., 2013). In order to explore the hydrodynamics of CFB risers, extensive computational fluid dynamics (CFD) studies have been carried out at different spatio-temporal scales (Van der Hoef et al., 2008; Ge et al., 2011), amongst of which continuum models are popular for industrial applications (Gidaspow, 1994; Jackson, 2000). At the same time, early studies have shown that particle clustering structures have a profound effect on the interphase drag force between gas and particles in gas–solid flow (Arastoopour and Gidaspow, 1979; Li et al., 1982, 1993; Dingrong and Yong, 1991; OBrien and Syamlal, 1993), which is a key input in continuum models. Therefore, it is necessary to consider the effect of particle clustering structures when modelling gas–solid flow. Various methods have been proposed to address this issue, such as the Energy Minimization Multi-Scale (EMMS) based approaches (Yang et al., 2003; Wang and Li, 2007; Wang et al., 2008, 2012; Zhou et al., 2014) and filtered two-fluid model (Agrawal et al., 2001; Igci et al., 2008; Parmentier et al., 2012; Milioli et al., 2013; Ozel et al., 2013; Schneiderbauer and Pirker, 2014; Zhu et al., 2016). The EMMS-based method has been applied to studying the heterogeneous gas–solid flow in CFB risers through two different ways: (i) modifying the interphase drag force in standard two-fluid model (Yang et al., 2003; Xiao et al., 2003; Wang and Li, 2007; Wang et al., 2008; Zhou and Wang, 2015; Zeneli et al., 2015; Liu et al., 2015); (ii) formulating a conceptually new two-fluid model (Wang et al., 2012), denoted as EMMS-based two-fluid model, using the principle of compromise of dominant mechanisms in competition originally proposed in the EMMS method (Li et al., 1988, 2013; Li and Kwauk, 1994). In the model, particle-rich dense phase and gas-rich dilute phase, corresponding to the physical realizations of particle-dominant mechanism and gas-dominant mechanism respectively, are treated as the interpenetrating continua (Wang et al., 2012), instead of treating gas and solid as the interpenetrating continua as in the standard two-fluid model (Gidaspow, 1994). Later on, it was shown that from practical viewpoint the model can be effectively simplified by assuming that all the particles are in the particle-rich dense phase (Zhou et al., 2014; Zhao et al., 2015). The EMMS-based two-fluid model is up to date established from the viewpoint of continuum mechanics, its microscopic foundation remains unknown. To this end, the statistical mechanics foundation of the simplified EMMS-based two-fluid model (Zhou et al., 2014; Zhao et al., 2015) is presented using generalized Boltzmann kinetic theory, which results in the two main contributions of present study: (i) A new method is developed by considering the fluctuations at different scales of gas phase simultaneously, with which we can for the first time derive the correct governing equations of gas phase via kinetic theory, in the sense that both the molecular stress and the Reynolds (or pseudoReynolds) stress can be obtained simultaneously, whereas all previous kinetic theory analyses failed to predict the Reynolds (or pseudo-Reynolds) stress due to the assumption of uniform structure, although this is physically always existent no matter how small the Reynolds number is; (ii) In case of particle phase, the generalized Boltzmann equation considering the spatio-temporal variation of the volume, density and velocity of clusters is firstly derived, a set of macroscopic transport equations are then derived in different phase spaces. It was shown that the governing equations of the dense phase in EMMS-based two-fluid model derived from continuum mechanics viewpoint corresponds to the macroscopic transport equations at (r, t ) space.

45

2. Kinetic theory analysis of gas phase Many studies have been devoted to the derivation of governing equations of gas–solid flow using Boltzmann kinetic theory for both gas and solid phases (Marble, 1963; Pai, 1977; Liu, 1987, 1993), unfortunately, the momentum conservation equation of gas phase obtained from all of those studies did not contain the Reynolds or pseudo-Reynolds stress term, which however should be appeared according to many other studies (Jackson, 2000; Ishii and Hibiki, 2011; Mehrabadi et al., 2015). The physical origin of pseudo-Reynolds stress is the presence of solid particles and/or particle clustering structures, which results in the variation of gas velocity around each particle or cluster, therefore, the pseudoReynolds stress is always existent no matter how small the Reynolds number is. In this section, we show how the governing equation of gas phase can be correctly derived from Boltzmann kinetic theory, by considering the fluctuations at molecular level and macroscopic level simultaneously. It turns out that the missing of pseudo-Reynolds stress term in previous studies is due to the neglect of macroscopic fluctuation. Note that previous studies (Marble, 1963; Pai, 1977; Liu, 1987, 1993) assume that the solid particles are homogeneously distributed in gas phase, whereas present study assumes that particle clusters are homogeneously distributed in gas phase, but this difference should have no influence on the appearance of pseudo-Reynolds stress term or not, although it does have an effect on the quantitative value of pseudo-Reynolds stress. Most of reported studies have modelled the turbulent flow in single phase and gas–solid two phase using macroscopic model, such as Reynolds Averaged Navier–Stokes (RANS) equations (Crowe et al., 1996; Pope, 2000). Recent studies on single phase turbulent flow have shown that there are two possible ways to consider the effect of turbulence in the framework of Boltzmann kinetic theory (Girimaji, 2007; Chen et al., 2003, 2004). In the works of Chen et al. (2003, 2004), the Boltzmann equation for velocity distribution function was filtered directly, the filtered Boltzmann equation was then used to model turbulent flow. However, according to the study of Girimaji (2007), there are three crucial difficulties associated with the studies of Chen et al. (2003, 2004): (i) the influence of the unresolved field on the filtered velocity distribution function evolution is implicitly considered via the averaged collision operator, that is, the influence of macroscopic turbulent fluctuation and the microscopic molecular collisions are both represented by the collision operator of the filtered Boltzmann equation. Therefore, it is difficult to explicitly derive the Reynolds stress term and then to utilize the abundant knowledge of RANS turbulent models, such as k − ϵ model; (ii) Since the collision operator should take both the microscopic molecular collision and macroscopic velocity fluctuation into consideration, it is not easy to find a suitable collision model; (iii) the final difficulty is related to find the proper filtered equilibrium distribution function. This represents an average over Maxwellian distribution function of different macroscopic velocities and cannot be Maxwellian about the averaged macroscopic velocity. To overcome the main difficulties related to the studies of Chen et al. (2003, 2004), Girimaji (2007) proposed a new method to separate the effect of macroscopic velocity fluctuation from molecular collision within the framework of Boltzmann kinetic theory via the so-called coordinate transformation technique. In this section, the method of Girimaji (2007) is extended to derive the correct governing equation of gas phase of gas–solid two-phase flow from Boltzmann kinetic theory. 2.1. Coordinate transform As has been discussed in the introduction, in the simplified

46

B. Zhao et al. / Chemical Engineering Science 156 (2016) 44–55

cell. The new gas molecule velocity becomes η = c − u . Because the new molecule velocity η is in one-to-one correspondence with the classical one c , the new probability density function g0(η, r, t|ugf ) = f0 (c , r, t|u′gf ), so the new velocity distribution can be expressed as:

g (η) ≡ (1 − f )g (η, r, t|ugf ) = (1 − f )f (c, r, t|u′gf ) = (1 − f )f (η + u , r, t|u′gf ).

(2)

2.2. The pseudo-Boltzmann equation and its filtering The deviation of the pseudo-Boltzmann equation for the evolution of the new velocity distribution function g (η ) of gas molecule is based on the classical Boltzmann equation. The relationship of spatio-temporal derivative between f (c) and g (η ) can be written as

Fig. 1. Schematic diagram of homogeneous distribution of clusters in a computational cell (the clusters are represented as black dots), which corresponds to the CFD cell used in the coarse grid simulation of the hydrodynamics of CFB risers. A point may be regarded as a cell used in direct numerical simulation of gas-cluster (or gas–particle) flow.

EMMS-based two-fluid model (Zhou et al., 2014), the dynamically heterogeneous gas–solid flow in CFB risers is simplified as clusters that are distributed uniformly in fluidizing gas, as shown in Fig. 1. Due to the disturbance of clusters on the gas flow, the gas velocities at different places are usually different, which is the source of the macroscopic velocity fluctuation and results in the appearance of Reynolds or pseudo-Reynolds stress term in the momentum equation, no matter how small the Reynolds number is. For example, the gas velocity near a cluster (point B) is different from the point away from the cluster (point A). Therefore, there are at least two intrinsic velocity fluctuations within the computational cell, one is the molecular scale velocity fluctuation, the other is the macroscopic velocity fluctuation. Note that previous studies (Marble, 1963; Pai, 1977; Liu, 1987, 1993) assumed that the macroscopic velocity at any points within a computational cell is not varying, meaning that there is no velocity fluctuation at macroscopic scale. In order to consider the effect of macroscopic velocity fluctuation and explicitly represent its effect, the coordinate transform technique of Girimaji (2007) is extended here. As in the kinetic theory of gases, we assume the velocity distribution function at any point (e.g. point A or point B) within a computational cell is Maxwellian (Chapman and Cowling, 1970; Kremer, 2010),

(

)

(

f (c) ≡ (1 − f )f c, r, t|u′gf = (1 − f )nf0 c, r, t|u′gf ⎛ ⎜ 1 = (1 − f )n⎜ ⎜⎜ 2π kT ⎝ m

)

(3)

∂g (η) ∂f (c) ∂g (η) ∂u = + · , ∂r ∂r ∂η ∂r

(4)

Following the classical kinetic theory, the standard Boltzmann equation for f (c) can be expressed as (Chapman and Cowling, 1970)

⎞ ∂f (c) ∂f (c) ∂ ⎛F + c· + ·⎜ f (c)⎟ = C1(f ) + C2(f ), ⎠ ∂t ∂r ∂c ⎝ m

(5)

where F denotes the external body force acting on the gas mole∂F cule, it is independent with the gas molecule velocity c , i.e., ∂c = 0, even if with the coordinate transformation, the expression of ∂F = 0 is also hold. Furthermore, C1(f ) represents the effect of the ∂η interaction between gas molecules, and the effect of the interaction between gas molecules and clusters is written as C2(f ) , which is the root of the interphase momentum exchange between the two phases. The substitution of (3) and (4) into the standard Boltzmann Eq. (5) leads to

⎡ ∂g (η) ⎡ ∂g (η) ∂g (η) ∂u ⎤ ∂g (η) ∂u ⎤ F ∂g (η) − · ⎥ + (η + u)·⎢ − · ⎥+ · ⎢ ⎣ ∂t ⎣ ∂r ∂η ∂t ⎦ ∂η ∂r ⎦ m ∂η = C1(g ) + C2(g ),

(6)

after rearranging, we have another form of (6) as

⎡ ∂g (η) ∂g (η) ⎤ ∂g (η) ∂g (η) + η· − a· = C1(g ) + C2(g ), ⎢ ⎥ + u· ⎣ ∂t ∂r ⎦ ∂r ∂η where a =

(

∂u ∂t

+

∂u η· ∂r

∂a ∂η

+

∂u u· ∂r

−

F m

∂u = ∂r . Furthermore, it is easy ∂g ∂(g u) ∂g ∂(g a) ∂u u· ∂r = ∂r − g ∂r and a· ∂η = ∂η

3

⎞2 ⎛ ⎞ ⎟ ⎜ (c − u′ )2 ⎟ gf ⎟ exp⎜− ⎟, kT ⎟⎟ ⎜⎜ ⎟⎟ 2 ⎠ ⎝ ⎠ m

∂g (η) ∂f (c) ∂g (η) ∂u = + · , ∂t ∂t ∂η ∂t

(7)

), so the partial derivation of a is to have the following equations: ∂a

− g ∂η =

∂(g a) ∂η

∂u

− g ∂r . Then, we have

another form of (6):

(1)

where u′gf and T are the macroscopic gas velocity and temperature at the point, respectively; k is the Boltzmann constant, m is the mass of a molecule and n is the number density of molecules; the occurrence of 1 − f (the volume fraction of gas phase) is due to the fact that there is only 1 − f probability is occupied by gas phase in two-phase flow. Let the reference frame moving with a velocity of u , which is equal to the macroscopic fluctuation velocity (u = u′gf − ugf ), where ugf is the macroscopic mean velocity of the computational

⎡ ∂g ∂g ⎤ ∂(g u) ∂(g a) + η· ⎥ + − = C1(g ) + C2(g ). ⎢ ⎣ ∂t ∂r ∂η ∂r ⎦

(8)

Eq. (8) is the pseudo-Boltzmann equation for describing the evolution of velocity distribution function g (η ). Note that we have dropping η in Eq. (8). Once obtaining the pseudo-Boltzmann equation for g (η ), it can be filtered and the filtered Boltzmann equation can be used to model the turbulence effectively (Girimaji, 2007; Chen et al., 2003, 2004). The filtered g (η ) is defined as g (η ), which is regarded as the velocity distribution function of gas of the computational cell. The

B. Zhao et al. / Chemical Engineering Science 156 (2016) 44–55

transport equation of g (η ) can be derived through the global average of (8), it is given as follows:

⎡ ∂g¯ ∂g¯ ⎤ ∂(g u ) ∂(g a ) + η· ⎥ + − = C1(g ) + C2(g ), ⎢ ⎣ ∂t ∂r ⎦ ∂r ∂η

(9)

which is the filtered pseudo-Boltzmann equation and is the key to the analysis of governing equations of gas phase via kinetic theory. It can be seen that with the coordinate transform technique, the macroscopic velocity fluctuation and the microscopic velocity fluctuation of molecules can be separated successfully. Specifically, the effect of macroscopic fluctuation is represented by the terms ∂(gu ) and ∂(ga) and the effect of microscopic fluctuation manifests ∂r

∂g¯

The transport equation of the zeroth moment and the first moment of the velocity η can be obtained by the integration of the filtered pseudo-Boltzmann Eq. (9) over the velocity domain. The mass and momentum densities of gas are respectively defined as

∫ mg (η)dη = (1 − f )mn = (1 − f )ρg ∫ mηg (η)dη = (1 − f )mnugf = (1 − f )ρg ugf ,

(10)

where the mass of the gas molecule is represented by m and the density of gas.

ρg is

2.3.1. Mass conservation equation Eq. (9) is multiplied by the molecule mass m and then the resulting equation is integrated over all the velocity η :

∫ m ∂∂gt dη + ∫ mη· ∂∂gr dη + ∫ m ∂(∂gru ) dη − ∫ m ∂(∂gηa ) dη ¯

a) dη = − (mg a ) ∫ m ∂(∂gηa ) dη = − ∫ ∂(mg ∂η

η =+∞

= 0, η =−∞

(15)

since it is assumed that there is no possibility that the velocity η is infinitely close to ±∞, so g ηη =+∞ = 0. =−∞ Furthermore, due to the conservation of mass during the collision between gas molecules,

∫ mC1(g )dη + ∫ mC2(g )dη

∂g¯ m dη = ∂t

∫ mC1(g )dη = 0,

(16)

it is also clear that the interactions between gas molecules and solid clusters have no influence on the total mass of the gas, so

∫ mC2(g )dη = 0.

(17)

Therefore, the mass conservation equation of gas can be written as

(18)

2.3.2. Momentum conservation equation In order to develop the momentum conservation equation, the Eq. (9) is multiplied by mη and then integrated, yielding

∫ mη ∂∂gt dη + ∫ mηη· ∂∂gr dη + ∫ mη· ∂(∂gru ) dη − ∫ mη· ∂(∂gηa ) dη ¯

¯

=

∫ mηC1(g )dη + ∫ mηC1(g )dη.

(19)

Via a similar procedure in the derivation of the mass conservation equation, the first and third term on the left hand side of Eq. (19) become

∫ mη ∂∂gt dη = ∂∂t [(1 − f )ρg ugf ], ¯

(20)

¯

∫ mη· ∂(∂gru ) dη = ∂∂r · ∫ (g u mη)dη + ∫ g u · ∂∂mrη dη (11)

In the following derivation, we note the independence in variables ( η , r , t), so the spatio-temporal derivative can be drawn out the integral of the gas molecule velocity η , and the mutual derivation of any two variables chosen from them equals zero identically. So the first and second terms on the left hand side of Eq. (11) are given by

∫

−

∂ ∂ [(1 − f )ρg ] + ·[(1 − f )ρg ugf ] = 0. ∂t ∂r

2.3. Governing equation

∫

where g (i) is a function of the gas molecule density ρg , the velocity ugf and temperature T, so it has no relation with the macroscopic fluctuation velocity u , i.e, g u = g¯ u = 0. The fourth term on the left hand side of Eq. (11) is

∂η

explicitly via the term η· ∂r and the source term C1(g ). It is interesting to note that as in the kinetic theory of gases or granular flows, at current stage it is not necessary to know the concrete form of the kernel for this filtration, because in this study we are only interested in the conservation equations from the generalized Boltzmann kinetic theory, however, different filtration kernels will result in different constitutive relationships, therefore, finding a suitable filtration kernel is the target of our future study.

=

47

∫

∂g¯ mη· dη = ∂r

∂(mg¯ ) ∂ dη = ∂t ∂t

∫

∂(mg¯ η) dη − ∂r

∫ ∫

∂ ¯ η = [(1 − f )ρg ], mgd ∂t ∂mη ∂ g¯ dη = ·[(1 − f )ρg ugf ]. ∂r ∂r

u) ∂ ∂ dη = · ∫ mg u dη = · ∫ mg¯ u¯ dη ∫ m ∂(∂gru ) dη = ∫ ∂(mg ∂r ∂r ∂r (14)

The reason for this derivation is that the distribution function of g (η ) can be expressed as an infinite series according the hypothesis of Chapman and Cowling (1970): g = g (0) + g (1) + g (2) + ⋯,

¯ ηη¯dη + ∫ g u ·0dη = 0. ∫ gm

(21)

∫ mηη· ∂∂gr dη = ∫ ∂∂r ·(mηηg¯ )dη − ∫ mg¯ ∂∂ηηr dη ¯

(12)

(13)

∂ · ∂r

The second term on the left hand side of Eq. (19) can be expressed as

=

From the definition of we have u′gf , u = u′gf − ugf = u′gf − ugf = ugf − ugf = 0 . Hence, the third term on the left hand side of Eq. (11) is

= 0.

=

∂ ⎡ ·⎢ ∂r ⎣ +2

¯ η + ∫ (η − ugf )(η − ugf )mgd ¯ η ∫ ugf ugf mgd ⎤

¯ η⎥ . ∫ ugf (η − ugf )mgd ⎦

(22)

The first term on the right hand side of Eq. (22) is

∂ · ∂r

⎛

⎞

∂ ¯ η= ¯ η · ugf ugf ∫ mgd ∫ ugf ugf mgd ⎠ ∂r ⎝ ⎜

=

⎟

⎤ ∂ ⎡ ·⎢ (1 − f )ρg ugf ugf ⎥. ⎦ ∂r ⎣

(23)

With the definition of the stress tensor σgf , the second term on the right hand side of Eq. (22) is

∂ · ∂r

∂ ¯ η=− ·[(1 − f )σgf ]. ∫ (η − ugf )(η − ugf )mgd ∂r

(24)

The third term on the right hand side of Eq. (22) can be expressed

48

B. Zhao et al. / Chemical Engineering Science 156 (2016) 44–55

as ∂ 2 · ∂r

∫

⎤ ∂ ⎡ ¯ η] = 2 ·⎢ u gf ηmgd ¯ η − u gf u gf mgd ¯ η⎥ u gf (η − u gf )mgd ⎦ ∂r ⎣ ⎤ ∂ ⎡ = 2 ·⎢ u gf (1 − f )ρg u gf − u gf u gf (1 − f )ρg ⎥ = 0. ⎦ (25) ∂r ⎣

∫

∫

Combining (23) and (24) with (25), we can obtain the final form of (22) as

⎡

⎤

∫ mηη· ∂∂gr dη = ∂∂r ·⎢⎣ (1 − f )ρg ugf ugf ⎥⎦ − ∂∂r ·[(1 − f )σgf ]. ¯

(26)

The fourth term on the left hand side of Eq. (19) can be written as

−

∫ mη· ∂(∂gηa ) dη = − ∫ ∂∂η ·(mηg a )dη + ∫ ga · ∂∂mηη dη =

⎛ ∂u F⎞ ∂u ∂u + mg ⎜ + η· + u· − ⎟ dη ⎝ ∂t m⎠ ∂r ∂r ∂u ∂u = 0 + mg dη + mgη· dη ∂t ∂r ∂u + mg u· dη − g F dη, ∂r

−

∫

∫

∫

∂u dη + ∂t

∫

mgη·

∂u dη = − ∂r

∫

mu

∂g dη − ∂t

∫

muη·

∂g dη . ∂r

(28)

The third term on the right hand side of Eq. (27) is

) ∂g u − ∫ mu · dη . ∫ mg u· ∂∂ur dη = ∫ m ∂(g∂uu r ∂r

(29)

So, the sum of the first three terms on the right hand side of Eq. (27) can be expressed as

∫ mg ∂∂ut dη + ∫ mgη· ∂∂ur dη + ∫ mg u· ∂∂ur dη ⎛

∂g ∂(g u) ⎟⎞ ∂(g uu) dη + m dη + ∂r ⎠ ∂r ∂r ⎛ ∂(g a) ⎞ ∂ mu ·⎜ + C1(g ) + C2(g )⎟dη + · mg uudη ⎝ ∂η ⎠ ∂r

=−

∫ mu·⎝ ∂∂gt

=−

∫

⎜

⎛ = − u·⎜ ⎝ ∂ + · ∂r

∫

⎞

a) dη + ∫ mC1(g )dη + ∫ mC2(g )dη⎟ ∫ ∂(mg ⎠ ∂η

∫ mg uudη

= − u ·{ (mg a) +

∫

+ η·

∂ ·(uu ∂r

η =+∞ η =−∞

+

∫ mC1(g )dη + ∫ mC2(g )dη}

∂ ¯ η) = ·[(1 − f )ρg uu]. ∫ mgd ∂r

(30)

Note that: (i) By the use of the pseudo-Boltzmann equation, it ∂g ∂g ∂(g u) ∂(g a) holds that ∂t + η· ∂r + ∂r = ∂η + C1(g ) + C2(g ). (ii) ∫ mC1(g )dη = 0 and ∫ mC2(g )dη = 0 have been established by the previous analysis about the mass conversation of gas. It can be seen that the pseudoReynolds stress [(1 − f )ρg uu] has been successfully obtained. The fourth term on the right hand side of Eq. (27) is

−

∫ g F dη = − (1 − f )ρg g,

(31)

Combination of (30) with (31) leads to the final form of (27):

−

∫ mη· ∂(∂gηa ) dη = ∂∂r ·[(1 − f )ρg uu] − (1 − f )ρg g.

∂ ·[(1 − f )ρg uu], ∂r

(34)

(27)

It is easy to show the first two terms on the right hand side of Eq. (27) lead to

mg

∂ ·[(1 − f )τgf ] ∂r

which consists of the pseudo-Reynolds stress term, the last term in the right hand side of Eq. (34), caused by macroscopic fluctuation of velocity due to the presence of particle clusters.

∫

∫

(33)

which represents the interphase momentum exchange rate between the gas phase and the cluster phase. Finally, with the decomposition of gas phase stress tensor σgf = − PI + τgf , we have the correct form of the momentum conservation equation

= − (1 − f )∇P + (1 − f )ρg g − Fdrag +

+∞ m(ηg a ) −∞

∫

∫ mηC2(g )dη = P∇(1 − f ) − Fdrag ,

⎤ ∂ ∂ ⎡ ((1 − f )ρg ugf ) + ·⎢ (1 − f )ρg ugf ugf ⎥ ⎦ ∂t ∂r ⎣

∫ ∂∂η ·(mηg a )dη + ∫ mgadη

=−

Because of the momentum conservation during the moleculemolecule interaction, the integral of the first source term ∫ mηC1(g )dη of the Eq. (19) equals to zero. And with the definition of C2(g ), the integral of the second source term can be expressed as (Enwald et al., 1996; Ishii and Mishima, 1984; Ishii, 1975),

(32)

3. Kinetic theory analysis of clustered dense phase The starting point of this section is the establishment of a generalized Boltzmann equation for the distribution function of clusters with size, density, velocity variation, which can reveal the spatio-temporal nature of the macroscopic structural characteristics from microscopic distribution function of clusters. As in the standard kinetic theory of gases, some macroscopic properties associated with the clustered dense phase are then defined, and transport equations of them at different phase spaces can be derived from the generalized Boltzmann equation. It can be shown that the macroscopic transport equations at (r, t ) space is in accordance with the governing equations of EMMS-based two-fluid model derived from continuum mechanics viewpoint (Zhou et al., 2014), thus laying a solid microscopic foundation of EMMS-based two-fluid model. Note that previous studies (Liu, 1987, 1993) have shown that the effect of the presence of clusters or particles on the gas phase can be implemented into Boltzmann equation in two different ways: One is that F only represents the effect of external body force (gravity in present study), then there is a term ∫ mηC2(g )dη representing the effect of the interaction between molecules and clusters (or particles), as has been shown in previous section. The other way is that the whole external forces (external body force plus interphase interaction force) can be represented by F as will be used in this section. 3.1. The generalized Boltzmann equation The first step in kinetic theory for clustered dense flow is the determination of how many independent variables should be used to specify the microscopic state of clusters. In the standard kinetic theory of gases, it usually assumes that with given translatory velocity c and position r of each gas molecule, the microscopic physical state of the studied system is known, therefore, the distribution function f (c , r, t ) only contains c , r , t as its independent variables (Chapman and Cowling, 1970; Kremer, 2010). However, for clusters in heterogeneous gas–solid flow, translatory velocity c and position r of each cluster are not sufficient to specify the microscopic state, even under the assumption of smooth and spherical clusters as in present study, because the size of cluster (or the volume of cluster) and the density of cluster (or the solid

B. Zhao et al. / Chemical Engineering Science 156 (2016) 44–55

concentration inside cluster) are varying, therefore, the volume of cluster v and the density of cluster ρ should also be added into the distribution function f as independent variables. That is, at given time t, a cluster can be specified by its volume v, density ρ, position r = (x1, x2, x3) and velocity c = (c1, c2, c3). Therefore, a cluster can be described by a point in an eight-dimensional space (v, ρ , c1, c2, c3, x1, x2, x3), and this space is called μ − space . Note that if we assume that the clusters are rough and non-spherical, the rotational motion of and the orientations of clusters should be included, which can in principle be done by adding more independent variables into the distribution function f. Following the standard process of kinetic theory, the generalized Boltzmann equation for the distribution function f can be derived as follows (Chapman and Cowling, 1970; Kremer, 2010): Supposing that f0 (v, ρ , c , r, t ) is the probability that a cluster is in the state (v, ρ , c , r) at the time t, the generalized distribution function f (v, ρ , c , r, t ) can then be defined as,

f (v, ρ , c, r, t ) = nf0(v, ρ , c, r, t ),

(35)

where, n ¼N/V denotes the number density of clusters in the system, N is the total number of clusters, V is total volume of this system. The number of clusters in the volume element dμ(t ) = dv dρ dcd r of μ − space at time t is

N (t ) = f (v, ρ , c, r, t )dv dρ dc dr.

N (t + Δt ) = f (v + Δv, ρ + Δρ , c + Δc, r + Δr, t + Δt )dμ(t + Δt ) (37)

.

Therefore, during the interval Δt , the change of the cluster number ΔN is,

ΔN = f (v + Δv, ρ + Δρ , c + Δc, r + Δr, t + Δt )dμ(t + Δt ) − f (v, ρ , c, r, t )dμ(t ).

(38)

The changes in the volume, density, velocity and position of the clusters during the time interval Δt are

Δv = GΔt ,

F Δc = Δt , m

Δr = cΔt ,

(39)

where F need to be stressed as the whole external forces imposed on the cluster, with

F m

=

dc ; dt

G is called by the growth rate and its

relationship with v can be expressed as G =

dv . dt

Furthermore, the

relationship between the two volume elements dμ(t + Δt ) and dμ(t ) is

dμ(t + Δt ) = |J|dμ(t ) = ∂(v(t + Δt ), ρ(t + Δt ), c(t + Δt ), r(t + Δt )) dμ(t ), ∂(v(t ), ρ(t ), c(t ), r(t ))

+

∂f ∂f ∂f ∂f Δv + Δρ + ·Δc + ·Δr ∂v ∂ρ ∂c ∂r

+

∂f Δt + O[(Δt )2]. ∂t

(42)

The substitution of (41) and (42) into (38) yields

⎤ ⎡ ∂f ΔN ∂ ∂ ⎛F ⎞ ∂ ⎛ dρ ⎞ ∂ =⎢ + ·(cf ) + ·⎜ f ⎟ + (Gf )⎥dμ(t ). ⎜ f⎟ + Δt ∂r ∂c ⎝ m ⎠ ∂ρ ⎝ dt ⎠ ∂v ⎦ ⎣ ∂t

(43)

On the other hand, r is used to describe the change of distribution function f due to the interaction between clusters, i.e. collisions. Hence, the change in number of cluster per unit of time is

ΔN = rdμ(t ). Δt

(44)

Therefore, we can deduce the final form of generalized Boltzmann equation as

∂f ∂ ∂ ⎛F ⎞ ∂ ⎛ dρ ⎞ ∂ + ·(cf ) + ·⎜ f ⎟ + (Gf ) = r , ⎜ f⎟ + ∂t ∂r ∂c ⎝ m ⎠ ∂ρ ⎝ dt ⎠ ∂v

(45)

Eq. (45) is a nonlinear integro-differential equation for the evolution of the volume, density and velocity distribution function of clusters. It is the central equation in kinetic theory analysis of the clustered solid phase.

⎛ dρ ⎞ F ∂⎜ ⎟ ∂ m ⎝ dt ⎠ ∂G ∂c |J | = 1 + Δt + Δt + Δt + Δt + O[(Δt )2]. ∂v ∂ρ ∂c ∂r

3.2. Definition of macroscopic variables The spatio-temporal development of heterogeneous clustering structure is the key to investigate hydrodynamics of fluidization in CFB risers. The detail derivations as follow will show that the macroscopic properties of cluster can be formulated by the generalized distribution function f (v, ρ , c , r, t ). First of all, two kinds of number densities of clusters at different phase spaces are defined as follows: the number density of clusters at (v, ρ , r, t ) space is

c (v, ρ , r, t ) = < f >c =

∫c fdc,

(46)

and the number density of clusters at (v, r, t ) space is

n(v, r, t ) = nv =

∫ρ < f >c dρ = ∫ρ ∫c fdcdρ.

(47)

Second, define ϕ(v, ρ , c , r, t ) as the property of the cluster with characteristic (v, ρ , c , r) at time t. By considering the velocity average of <ϕ >c (v, ρ , r, t ) = < ϕ>c , it follows that

<ϕ >c (v, ρ , r, t ) = < ϕ >c =

∫c ϕfdc ∫c fdc

, (48)

which represents the average value of ϕ at (v, ρ , r, t ) space. Then conducting the velocity and density average at the same time, which leads to

<ϕ>(v, r, t ) = < ϕ > =

(40)

∫ρ < ϕ >c < f >c dρ ∫ρ < f >c dρ

=

∫ρ ∫c ϕfdcdρ ∫ρ ∫c fdcdρ

. (49)

Next let ϕ be specific, the average velocity u v(v, r, t ) = u v of the clusters with the volume v is obtained by the substitution of ϕ = c into (49)

and the Jacobian matrix is

( )

Furthermore, expanding in Taylor series, we have

f (v + Δv, ρ + Δρ , c + Δc, r + Δr, t + Δt ) = f (v, ρ , c, r, t )

(36)

Similarly, the number of clusters in the volume element dμ(t + Δt ) of μ − space at time t + Δt is

dρ Δρ = Δt , dt

49

(41)

f (v + Δv, ρ + Δρ , c + Δc , r + Δr, t + Δt )

u v (v, r, t ) = u v = ⟨c⟩ =

∫ρ < c >c < f >c dρ ∫ρ < f >c dρ

=

∫ρ ∫c cfdcdρ ∫ρ ∫c fdcdρ

. (50)

Similarly, the average density ρv (v, r, t ), mass 〈mv 〉(v, r, t ) and

50

B. Zhao et al. / Chemical Engineering Science 156 (2016) 44–55

momentum 〈mv c〉 of the clusters with the volume v are defined as follows:

ρv (v, r, t ) = ρv = < ρ > =

∫ρ < ρ >c < f >c dρ ∫ρ < f >c dρ

∫ρ ∫c mvfdcdρ

〈mv 〉(v, r, t ) = 〈mv 〉 =

=

∫ρ ∫c fdcdρ

=

∫ρ ∫c ρfdcdρ ∫ρ ∫c fdcdρ

v ∫ ∫ ρfdcdρ ρ c

∫ρ ∫c fdcdρ

, (51)

∫c ϕ ∂∂ft dc = ∫c ∂∂t (ϕf )dc − ∫c f ∂∂ϕt dc = ∂∂t ∫c ϕfdc − ∫c ∂∂ϕt fdc

= vρv , (52)

〈mv c〉(v, r, t ) = 〈mv 〉〈c〉 = vρv u v .

(53)

Third, at position r and at time t, the macroscopic properties of the clustered solid phase is given by studying the ergodic of all features of clusters one by one. With the combination of the definition of (47), (50)–(53), the volume fraction of the clustered solid phase is given by

f (r, t ) = f =

where the suffix c on the integral indicates that it's taken over the whole of velocity domain. Because of the independence of t and c , the time derivatives can be drawn out the integral and with (49) about the definition of <•>c , the first term on the left hand side of Eq. (60) is

∫v vnvdv.

=

∂ ∂ (c < ϕ>c) − < f >c < ∂ϕt >c . ∂t

(61)

Moreover, other terms on the left hand side of Eq. (60) can be deduced as follows:

∫c ϕ ∂∂r ·(cf )dc = ∫c ∂∂r ·(fϕc)dc − ∫c f c· ∂∂ϕr dc ∂ ∂ϕ · fϕcdc − f c· d c c ∂r c ∂r ∂ ∂ϕ = ·(c < ϕc>c ) − < f >c < c· >c , ∂r ∂r

∫

=

(54)

∫

(62)

The mass of the clustered solid phase per unit of system volume is

f (r, t )ρdense (r, t ) = fρdense =

⎛

∫v 〈mv〉nvdv = ∫v vρv nvdv.

(55)

⎞

The momentum of the clustered solid phase per unit of system volume is

∫v 〈mvc〉nvdv = ∫v vρv uvnvdv.

∫v 〈mvc〉nvdv ∫v 〈mv〉nvdv

=

∫v vρv uvnvdv ∫v vρv nvdv

∫c d⎜⎝ mF fϕ⎟⎠ − < f >c < mF · ∂∂ϕc

(58)

and at (v, r, t ) space it is

⎞

(59)

3.3. Derivation of macroscopic transport equation The integration of the generalized Boltzmann Eq. (45) over different domains (velocity, density and volume domains) one by one will lead to the transport equations for the property ϕ(v, ρ , c , r, t ) of cluster in different spaces (Chapman and Cowling, 1970; Kremer, 2010).

=

⎛

⎞

⎛

+

∫c

∫c ϕrdc,

F ∂ϕ · >c m ∂c

− < f >c < −∞

F ∂ϕ · >c , m ∂c

⎛

⎞

⎛

⎞

(63) ⎛

⎞

∫c ∂∂ρ ⎜⎝ϕf ddtρ ⎟⎠dc − < f >c < ddtρ ∂∂ϕρ

>c ,

(64)

∂ ∂ϕ Gf dc (Gfϕ)dc − c ∂v c ∂v ∂ ∂ϕ = < f >c < Gϕ >c − < f >c < G >c . ∂v ∂v

∫

∫

(65)

With the comprehensive treatment of all terms above and rearranging appropriately, we finally obtain the transport equation of <ϕ>c at (v, ρ , r, t ) space ∂ (⟨f ⟩c ⟨ϕ⟩c ) + ∂∂r ·(⟨f ⟩c ⟨ϕc⟩c ) + ∂t ∂ϕ ∂t

+ ⟨f ⟩c

+ ⟨f ⟩c c· c

dρ ∂ϕ dt ∂ρ

⎛

∂ϕ ∂r

+ ⟨f ⟩c c

⎞

∫c ∂∂ρ ⎜⎝ϕf ddtρ ⎟⎠dc + ∂∂v + ⟨f ⟩c c

∂ϕ G ∂v

+ c

< f >c < Gϕ >c

F ∂ϕ · m ∂c

∫c ϕrdc.

c

(66)

⎞

∫c ϕ ∂∂ft dc + ∫c ϕ ∂∂r ·(cf )dc + ∫c ϕ ∂∂c ·⎜⎝ mF f ⎟⎠dc + ∫c ϕ ∂∂ρ ⎜⎝ ddtρ f ⎟⎠dc ∂ ϕ (Gf )dc = ∂v

>c

∫c ϕ ∂∂v (Gf )dc = ∫c ∂∂v (Gfϕ)dc − ∫c Gf ∂∂ϕv dc

= ⟨f ⟩c

3.3.1. Transport equation at (v, ρ , r, t ) space To get the transport equation of <ϕ >c (v, ρ , r, t ) at (v, ρ , r, t ) space, the generalized Boltzmann Eq. (45) need to be multiplied by ϕ(v, ρ , c , r, t ), and then integrated over the velocity domain, i.e.,

>c

∫c ϕ ∂∂ρ ⎜⎝ ddtρ f ⎟⎠dc = ∫c ∂∂ρ ⎜⎝ϕf ddtρ ⎟⎠dc − ∫c f ⎜⎝ ddtρ ∂∂ϕρ ⎟⎠dc

=

wv(v, r, t ) = wv = 〈C〉 = u v − udense .

+∞

= − < f >c <

What is more, according to (50) and (57), the definition of two fluctuation velocities of the cluster at different spaces are given below: at (v, ρ , c , r, t ) space it is

C(v, ρ , c, r, t ) = C = c − udense ,

⎞

⎛

⎛F ⎞ = ⎜ f ϕ⎟ ⎝m ⎠

⎛

(57)

⎞

=

(56)

.

⎛

∫c ∂∂c ·⎜⎝ mF fϕ⎟⎠dc − < f >c < mF · ∂∂ϕc

The velocity of the clustered solid phase is computed by combining (55) with (56), i.e.,

udense (r, t ) = udense =

⎞

=

f (r, t )ρdense (r, t )udense (r, t ) = fρdense udense =

⎛

∫c ϕ ∂∂c ·⎜⎝ mF f ⎟⎠dc = ∫c ∂∂c ·⎜⎝ mF fϕ⎟⎠dc − ∫c f mF · ∂∂ϕc dc

(60)

3.3.2. Transport equation at (v, r, t ) space The transport equation for 〈ϕ〉 at (v, r, t ) space can be obtained by integrating both sides of (66) over density domain

B. Zhao et al. / Chemical Engineering Science 156 (2016) 44–55

⎛

⎞

integral of (71) over the cluster volume domain (0, + ∞). Let ϕ = mv and integrate both side of (71), we have

∫ρ ∂∂t (c ⟨ϕ⟩c)dρ + ∫ρ ∂∂r ·(⟨f ⟩c⟨ϕc⟩c)dρ + ∫ρ ∫c ∂∂ρ ⎜⎝ϕf ddtρ ⎟⎠dcdρ +

∫ρ ∂∂v ⟨f ⟩c⟨Gϕ⟩cdρ

∫

∫

∫

∫

∫

∫∫

(67)

∫ρ ∂∂t ( < f >c < ϕ >c )dρ = ∂∂t ∫ρ ( < f >c < ϕ >c ) ∂ (nv < ϕ>). ∂t

(68)

Since the cluster with density ρ =+∞ does not appear in the system, i.e., f (v, + ∞ , c , r, t ) = 0. Besides, the fact that the cluster must be composed by particles leads to ρ ≠ 0, i.e., f (v, 0, c , r, t ) = 0. So, the third term on the left hand side of the Eq. (67) can be expressed as

⎛

⎞

⎡

⎝

⎠

⎣

⎛

⎞

⎤

∫ρ ∫c ∂∂ρ ⎜⎜ ϕf ddtρ ⎟⎟dcdρ = ∫c ⎢⎢ ∫ρ ∂∂ρ ⎜⎝ ϕf ddtρ ⎟⎠dρ⎥⎥dc =

∫c [(ϕf ddtρ )

⎦

+∞

]dc = 0. 0

(69)

It's worth stressing again r describes the change of f due to the collision with clusters in the generalized Boltzmann equation, and ϕrdcdρ measures the change of ϕ due to the change of cluster number by the collision. The double integral of ϕrdcdρ at cluster velocity space and density space represents the change of the common property ϕ of all cluster with the volume v. The sixth term on the right hand side of the Eq. (67) thus becomes

∫ρ ∫c ϕrdcdρ = Γ〈nvϕ〉,

(70)

where Γ is change operator. Dealing with other terms in Eq. (67) by the similar method of (68), the final form of (67) is given by

∫

∫

∫

∫

∫

∫

which is the continuity equation for the clustered solid phase. To simplify the notation, it is necessary to do the derivation further. During the derivation, in addition to using the definitions of f, udense , fρdense and fρdense udense , we need to note the mutual derivation of any two variables chosen in v, ρ , c , r, t equals zero identically by the microscopic independence, and with the same reason, the derivatives and the integral in volume domain can exchange operation sequence. Hence, (73) can be rewritten as

∂ ∂ ·(fρ )+ ·(fρ udense ) = 0. ∂t dense ∂r dense It should be note

∞ ∂ (nv〈mv G〉)dv 0 ∂v

∫

(74) = (nvG〈mv 〉)

+∞ 0

does not appear,

since the volume v of cluster must no less than the volume of single particle, i.e, n0 = 0, and it seems reasonable to assume that it is impossible the volume v of cluster is infinitely close to +∞, i.e, n+∞ = 0. Besides, by considering the source term caused by collision,

∞

∫0 Γ〈mvnv〉dv = Γ (fρdense ) is 0 due to the conservation of mass

of clusters during cluster-cluster interactions. Note that in the formulation of continuity equation for aggregating particles in the work of Liu (2011), a similar treatments were used, whose approach was multiplying mv to both sides of (72) and then integrating over the volume domain. This method must base on the assumption that the mass of particles are constant. While, in fact, once aggregation and breakage of particles taking place, the mass of particles should be changed. To get the momentum equation of clustered solid phase, the substitution of ϕ = mv c into (71) and the integration of it leads to ∞ ∂ ∂ (nv 〈mv c〉)dv + ·(nv 〈mv cc〉)dv 0 ∂r ∂t ∞ ∂ + (nv 〈Gmv c〉)dv 0 ∂v ∞ ∞ ∂(mv c) ∂(mv c) = nv〈 〉dv + nv 〈c· 〉dv 0 0 ∂t ∂r ∞ ∞ F ∂(mv c) dρ ∂(mv c) nv〈 nv 〈 + · 〉dv + 〉dv 0 0 dt ∂ρ mv ∂c ∞ ∞ ∂(mv c) nv 〈G Γ 〈nv mv c〉dv. + 〉dv + 0 0 ∂v

∫0

∞

∫

∫

(71)

Noting that setting ϕ = 1 in (71), the standard population balance equation can be obtained:

∂ ∂ ∂ (nv ) + ·(nv 〈c〉) + (nv 〈G〉) = Γ 〈nv 〉. ∂t ∂r ∂v

∫

∫

∂ ∂ ∂ ∂ϕ (nv < ϕ>) + ·(nv ⟨ϕc⟩) + (nv ⟨Gϕ⟩) = nv ⟨ ⟩ ∂t ∂r ∂v ∂t ∂ϕ F ∂ϕ ∂ϕ dρ ∂ϕ ⟩ + nv ⟨c· ⟩ + nv ⟨ · ⟩ + nv ⟨G ⟩ + nv ⟨ ∂r ∂v m ∂c dt ∂ρ + Γ ⟨nv ϕ⟩.

∞

∫

where the suffix ρ on the integral indicates that it's taken over the whole of density domain. Following the procedure of (61), the first term on the left hand side of the Eq. (67) is given by

dρ =

∞ ∂ ∞ ∂ ∂ (nv 〈mv 〉)dv + ·(nv 〈mv c〉)dv + (nv 〈Gmv 〉)dv 0 ∂r 0 ∂v ∂t ∞ ∞ ∞ ∂m F ∂mv dρ ∂mv = nv 〈G v 〉dv + nv 〈 · 〉dv + nv 〈 〉dv 0 0 0 ∂v m ∂c dt ∂ρ ∞ ∞ ∞ ∂mv ∂m nv 〈c· nv 〈 v 〉dv + Γ 〈nv mv 〉dv, + 〉dv + (73) 0 0 0 ∂r ∂t

∫0

F ∂ϕ ∂ϕ ∂ϕ = ⟨f ⟩c ⟨ ⟩c dρ + ⟨f ⟩c ⟨c· ⟩c dρ + ⟨f ⟩c ⟨ · ⟩c dρ ρ ρ ρ m ∂c ∂t ∂r dρ ∂ϕ ∂ϕ ϕrdcdρ , + ⟨f ⟩c ⟨ ⟩c dρ + ⟨f ⟩c ⟨G ⟩c dρ + ρ ρ ρ c dt ∂ρ ∂v

51

(72)

Present study focuses on the derivation of conservation equations from generalized Boltzmann kinetic theory, therefore, the explicit from of r is not needed. However, in the future a suitable physical model for the interaction between clusters (i.e. a model for cluster aggregation and breakage) will be established, Eq. (72) then can be used to predict the size distribution of clusters in gas–solid flow. With the use of the ratio of the third moment to secondary moment, the prediction of the spatio-temporal variation of cluster size can be achieved, which is very helpful to revise the gas–solid drag force. 3.3.3. Transport equation at (r, t ) space Similar to the deviation of the transport Eq. (71) at (v, r, t ) space, the transport equation at (r, t ) space is deduced by the

∫

∫

∫

∫

∫

(75)

The second and third terms on the left hand side of Eq. (75) deserve special attention:

∫0

∞

⎞ ∂ ∂ ⎛⎜ ∞ ·(nv 〈mv cc〉)dv = · nv 〈mv cc〉dv⎟ ⎠ ∂r ∂r ⎝ 0 ∂ = ∂r ⎛ ∞ ⎞ ·⎜ nv 〈mv (C + udense )(C + udense )〉dv⎟ ⎝ 0 ⎠

∫

∫

=

⎞ ∂ ⎛⎜ ∞ · nv 〈mv CC〉dv⎟ ⎠ ∂r ⎝ 0 ⎞ ∂ ⎛⎜ ∞ + · nv 〈mv udense udense 〉dv⎟ ⎠ ∂r ⎝ 0 ⎞ ∂ ⎛⎜ ∞ +2 · nv 〈mv Cudense 〉dv⎟. ⎠ ∂r ⎝ 0

∫

∫

∫

The first term on the right hand side of Eq. (76) is

(76)

52

B. Zhao et al. / Chemical Engineering Science 156 (2016) 44–55

∂ ⎛⎜ · ∂r ⎝

∫0

∞

⎞ ∂ ·[ − pdense I¯ + fτdense ] nv 〈mv CC〉dv⎟ = − ⎠ ∂r ∂ =− ·(fσdense ). ∂r

(77)

The second term on the right hand side of Eq. (76) is

∂ ⎛⎜ · ∂r ⎝

∫0

∞

∞ ⎞ ∂ ·(u nv 〈mv udense udense 〉dv⎟ = u nv 〈mv 〉dv) ⎠ ∂r dense dense 0 ∂ = ·(fρ udense udense ). ∂r dense

∫

(78)

Combined with the definition of the velocity (57) and the fluctuation velocity (58) of the clustered solid phase, the third term on the right hand side of Eq. (76) can be written as 2

∂ ⎛⎜ · ∂r ⎝

∞ ⎞ ⎞ ∂ ⎛ n v〈mvCu dense〉dv⎟ = 2 ·⎜ u dense n v〈mv(c − u dense)〉dv⎟ ⎠ ⎠ 0 ∂r ⎝ ∞ ∞ ⎞ ∂ ⎛⎜ = 2 · u dense n v〈mvc〉dv − u dense n v〈mvu dense〉dv⎟ ⎠ 0 0 ∂r ⎝ ∂ = 2 ·(u dense(fρdense u dense) − u denseu densefρdense ) = 0. ∂r

∫0

∞

∫

∫

∫

(79)

The third term on the right hand side of Eq. (75) is given by

∫0

∞

nv 〈

F ∂(mv c v) · 〉dv = mv ∂c v

∫0

∞

nv〈

F ·mv 〉dv = mv

∫0

∞

nv 〈F〉dv,

(80)

where 〈F〉 represents the average external force imposed on the clusters with the volume v. In general, it mainly includes the pressure gradient force −v∇P , the gravity 〈mv 〉g and the drag force 〈Fdrag 〉, and other forces (such as the added mass force, the Basset force, the Magnus effect and the Saffman lift force) can be negligible. Thus we have

∫0

∞

nv 〈F〉dv =

∫0

∞

nv ( − v∇P + 〈mv 〉g + 〈Fdrag 〉)dv

= − f ∇P + fρdense g + Fdrag .

(81)

Therefore, the final form of (75) can be obtained as follows:

∂ ∂ ·(fρ ·(fρ udense ) + udense udense ) ∂t dense ∂r dense ∂ = − f ∇P + ·(fσdense ) + fρdense g + Fdrag . ∂r

(82)

Note that: (i) There is a slight difference between (82) and our previous work (Zhou et al., 2014) in term of the stresses of clus∂ tered solid phase, owing to the fact that the effect of ∂r ·( − pdense I¯) on the simulation results is minor (Gidaspow, 1986), which has been ignored in our previous work (Wang et al., 2012; Zhou et al., 2014). (ii) Because of the conservation of the total momentum of ∞ clusters, the source term ∫ Γ 〈nvmv c〉dv = Γ (fρdense udense) ≡ 0 whe0 ther collisions occurring or not; (iii) The simplification of other items in (75) is analogous to the previous analysis of the continuity equation; and (iv) in our previous studies, the bulk viscosity of the clustered solid phase is set to be zero due to the usage of empirical correlation, but now it is possible to derive the shear and bulk viscosity from first principle.

4. CFD simulation As has been indicated that one of main contribution is the derivation of pseudo-turbulent term from kinetic theory, which has been neglected in previous studies (Zhou et al., 2014; Zhao et al., 2015). Therefore, we choose the same riser (shown in Fig. 2) discretized using non-uniform hexahedral grid (dimensionless sizes: 47∼254) with the same operational conditions (Pärssinen and

Fig. 2. The geometry of simulated risers (for a clear show, the riser is not at scale).

Zhu, 2001a, 2001b) as in our previous works to investigate the influence of gas-phase pseudo-turbulence, where the corresponding parameters used in numerical simulations are summarized in Table 1. Note that the simulation results obtained from the EMMS-based two-fluid model without the effect of gas-phase pseudo-turbulence and standard two-fluid model have been published in our previous studies (Zhou et al., 2014; Zhao et al., 2015). By virtue of the gas-phase pseudo-Reynolds stress included in the momentum conservation Eq. (34), it is necessary to select a reasonable constitutive relation to close it. Here we employ the analogous Bossinesq hypothesis to relate the pseudo-Reynolds stress (B.8) to the macroscopic gas-phase velocity gradient. Therefore, we can close the pseudo-Reynolds stress by defining the effective viscosity μ⁎ (B.5) as the summation of the true air viscosity μg and the pseudo-turbulent viscosity μpt (B.4), and its form is assumed as the same with the one in a suspension of forcefree spherical particles in Newton ambient fluid (Batchelor and Green, 1972). All the constitutive relations, including the constitutive relations of cluster phase and the interphase drag coefficient, are summarized in Table 2, which are still obtained from empirical correlations as previous studies (Zhou et al., 2014; Zhao et al., 2015). Figs. 3 and 4 indicate that the inclusion of the gas phase pseudo-turbulence has almost no influence on the axial distributions of solids concentration and the solids inventory in the whole bed, whereas the small existing differences is most likely due to the relatively short statistical time for calculating mean values (10 s). It can also be seen that the standard two-fluid model, which stands for the two-fluid model available in commercial software with kinetic theory of granular flow for particulate phase stress and the Gidaspow's drag correlation for interphase drag force, significantly low-predicts the solids inventory and axial solids concentration. Whereas the EMMS-based two-fluid models

B. Zhao et al. / Chemical Engineering Science 156 (2016) 44–55

53

Table 1 Physical properties and parameters used in simulations. 1500

Particle density, kg/m3 Particle diameter, μm

67 1.2

Air density, kg/m3 Air viscosity, Pa s

1.78 × 10−5 8

Gas velocity, m/s

300

Solid flux, kg/(m2 s) Grid number

5.5 550

300

195,000

Table 2 Summary of constitutive relations for the EMMS-based two fluid model. Mixture velocity of dense phase:

u dense =

εgc ρg u gc + (1 − εgc )ρs u sc εgc ρg + (1 − εgc )ρs

(B.1)

Fig. 3. Comparison of the volume-averaged solids concentration of bed simulated by EMMS-based two-fluid model and standard two-fluid model with experimental data of Pärssinen and Zhu (2001a).

Densityof dense phase:

ρdense = εgc ρg + (1 − εgc )ρs

(B.2)

Viscosity of dense phase: −2.5εs , max ⎛ εs, dense ⎞ ⎟⎟ μdense = μg ⎜⎜ 1 − εs, max ⎠ ⎝

with

εs, max = 0.63

(B.3)

Pseudo- turbulent viscosity of gas phase: 2

μ pt = (2.5f + 7.6f )μg

(B.4)

Effective viscosity of gas phase: 2

μ⁎ = μg + μ pt = (1 + 2.5f + 7.6f )μg

(B.5)

Stressof each phase: T )− τdense = μdense (∇u dense + ∇u dense

τgf = μg (∇u gf + ∇uTgf ) −

2 μ (∇u dense)I 3 dense

2 μ (∇u gf )I 3 g

(B.6)

(B.7)

Pseudo- Reynolds stress of gas phase:

−ρg uu = μ pt (∇u gf + ∇uTgf ) −

2 μ (∇u gf )I 3 pt

(B.8)

Drag coefficient:

⎛ ⎞ ρ f (1 − f )|u gf − u dense| g 17.3 β = ⎜⎜ + 0.336⎟⎟ (1 − f )−2.8 Re dcl ⎝ ⎠

Re =

(B.9)

5. Conclusion

(1 − f )ρg |u gf − u dense|dcl μg

uniform velocity profile at the inlet. Fig. 5, for example, compares the simulated radial profiles of solids concentrations with experimental data at various heights. It can be seen that the basic core–annulus structure has been captured, and the solid concentration at the centre can be predicted a litter bit more accurately using the EMMS-based two-fluid model with gas phase pseudo-turbulence. But as a whole, the solids concentration near the wall is low-predicted and the solids concentration at the centre is over-predicted, meaning that EMMSbased two-fluid model predicts a flatter radial solids concentration profile. Those should be other indications of the requirement of improving the prediction of cluster size, since previous study has shown that cluster size is the most important input parameter in the EMMS-based two-fluid model. As have been indicated, a work based on present study is ongoing to formulate a suitable model for cluster breakup and aggregation and then to predict a more reasonable cluster size, which are expected to improve the results of CFD simulations. Fig. 6 shows the effect of gas phase pseudo-turbulence on the axial solid velocity profiles. It can be seen that (i) the particle velocity profiles can be predicted reasonably well, including the basic feature of dense suspension upflow regime, that is, the particle velocities at the near wall region are upflow, instead of downflow as in the classical CFB risers; (ii) gas phase pseudo-turbulence has a notable effect of the axial solid velocity profiles, especially at the center region where the solids concentration is relatively low. In conclusion, these two aspects mean that gas phase pseudo-turbulence is important in cases of dilute gas–solid flow.

(B.10)

predict a much better agreement with experimental data, although the solids concentration at the upper part of the riser is over-predicted and the solids concentration near the inlet is lowpredicted. The former deviation is most possibly an indication that the cluster size in the upper part is over-predicted by EMMS model, and the latter deviation is due to the assumption of

A statistical mechanics foundation of EMMS-based two-fluid model has been established. It features the usage of filtered pseudo-Boltzmann equation and the coordinate transform technique in the derivation of governing equations of gas phase, and a generalized Boltzmann equation in the derivation of governing equations of cluster phase. CFD simulations with empirical correlations for the effective viscosity of both gas and cluster phases and EMMS-predicted cluster size for the interphase drag force indicate that EMMS-based two-fluid model can predict the hydrodynamics of high-density risers reasonably. The key deficiency of EMMS-based two-fluid model has also been highlighted, that is

54

B. Zhao et al. / Chemical Engineering Science 156 (2016) 44–55

Fig. 4. Comparison of predicted axial solids concentration profiles by EMMS-based two-fluid model with/without the effect of gas-phase pseudo-turbulence and standard two-fluid model with experimental data of Pärssinen and Zhu (2001a).

Fig. 5. Comparison between numerically simulated by EMMS-based two-fluid model and experimentally measured radial solid concentration profiles at different heights of the bed.

the requirement of a better prediction of spatio-temporal variation of cluster size. Therefore, our next step is to establish physical models for the breakage and aggregation of clusters. With the established breakage and aggregation kernels, the spatio-temporal variation of cluster size and then a better interphase drag force can be determined from the population balance equations derived in this study. Furthermore, as in the kinetic theory of gases and granular flows, the particulate phase stresses (such as shear and viscosities) can in principle be predicted from the generalized Boltzmann equation using Chapman–Enskog expansion. All those models together will result in a self-consistent and complete formation of EMMS-based two-fluid model for heterogeneous gas– solid flow in CFB risers with dynamically spatio-temporal variation of clustering characteristics.

Fig. 6. Comparison of axial particle velocity between simulation and experiment at different heights of the bed.

Acknowledgements This study is financially supported by the National Natural Science Foundation of China under the Grant Nos. 21206170 and 21422608, and the “Strategic Priority Research Program” of the Chinese Academy of Sciences under the grant No. XDA07080200.

B. Zhao et al. / Chemical Engineering Science 156 (2016) 44–55

References Agrawal, K., Loezos, P.N., Syamlal, M., Sundaresan, S., 2001. The role of meso-scale structures in rapid gas–solid flows. J. Fluid Mech. 445, 151–185. Arastoopour, H., Gidaspow, D., 1979. Analysis of igt pneumatic conveying data and fast fluidization using a thermohydrodynamic model. Powder Technol. 22 (1), 77–87. Batchelor, G., Green, J., 1972. The determination of the bulk stress in a suspension of spherical particles to order c 2. J. Fluid Mech. 56 (03), 401–427. Berruti, F., Pugsley, T., Godfroy, L., Chaouki, J., Patience, G., 1995. Hydrodynamics of circulating fluidized bed risers: a review. Can. J. Chem. Eng. 73 (5), 579–602. Chapman, S., Cowling, T.G., 1970. The Mathematical Theory of Non-Uniform Gases. Cambridge University Press, New York. Chen, H., Kandasamy, S., Orszag, S., Shock, R., Succi, S., Yakhot, V., 2003. Extended Boltzmann kinetic equation for turbulent flows. Science 301 (5633), 633–636. Chen, H., Orszag, S.A., Staroselsky, I., Succi, S., 2004. Expanded analogy between Boltzmann kinetic theory of fluids and turbulence. J. Fluid Mech. 519, 301–314. Cocco, R., Shaffer, F., Hays, R., Karri, S.R., Knowlton, T., 2010. Particle clusters in and above fluidized beds. Powder Technol. 203 (1), 3–11. Crowe, C., Troutt, T., Chung, J., 1996. Numerical models for two-phase turbulent flows. Annu. Rev. Fluid Mech. 28 (1), 11–43. Dingrong, B., Yong, J., 1991. The interaction between gas and particles in a vertical gas–solid flow system. J. Chem. Ind. Eng. (China) 42 (6), 697–703 (In Chinese). Enwald, H., Peirano, E., Almstedt, A.-E., 1996. Eulerian two-phase flow theory applied to fluidization. Int. J. Multiph. Flow 22, 21–66. Ge, W., Wang, W., Yang, N., Li, J., Kwauk, M., Chen, F., Chen, J., Fang, X., Guo, L., He, X., et al., 2011. Meso-scale oriented simulation towards virtual process engineering (vpe) the emms paradigm. Chem. Eng. Sci. 66 (19), 4426–4458. Gidaspow, D., 1986. Hydrodynamics of fiuidizatlon and heat transfer: supercomputer modeling. Appl. Mech. Rev. 39 (1), 1–23. Gidaspow, D., 1994. Multiphase flow and fluidization: Continuum and Kinetic Theory Descriptions. Academic Press, Boston. Girimaji, S.S., 2007. Boltzmann kinetic equation for filtered fluid turbulence. Phys. Rev. Lett. 99 (3), 034501. Harris, A., Davidson, J., Thorpe, R., 2002. The prediction of particle cluster properties in the near wall region of a vertical riser (200157). Powder Technol. 127 (2), 128–143. Igci, Y., Andrews, A.T., Sundaresan, S., Pannala, S., O'Brien, T., 2008. Filtered twofluid models for fluidized gas-particle suspensions. AIChE J. 54 (6), 1431–1448. Ishii, M., 1975. Thermo-fluid Dynamic Theory of Two-phase Flow. NASA STI/Recon Technical Report A 75, 29657. Ishii, M., Hibiki, T., 2011. Thermo-Fluid Dynamics of Two-Phase Flow, Second edition. Springer, New York. Ishii, M., Mishima, K., 1984. Two-fluid model and hydrodynamic constitutive relations. Nucl. Eng. Des. 82 (2), 107–126. Jackson, R., 2000. The Dynamics of Fluidized Particles. Cambridge University Press, Cambridge. Kremer, G.M., 2010. An Introduction to the Boltzmann Equation and Transport Processes in Gases. Springer Science & Business Media. Kunii, D., Levenspiel, O., 1991. Fluidization Engineering. Butterworth-Heinemann, Boston. Lackermeier, U., Rudnick, C., Werther, J., Bredebusch, A., Burkhardt, H., 2001. Visualization of flow structures inside a circulating fluidized bed by means of laser sheet and image processing. Powder Technol. 114 (1), 71–83. Li, H., Xia, Y., Tung, Y., Kwauk, M., 1991. Micro-visualization of clusters in a fast fluidized bed. Powder Technol. 66 (3), 231–235. Li, H., Zhu, Q., Liu, H., Zhou, Y., 1995. The cluster size distribution and motion behavior in a fast fluidized bed. Powder Technol. 84 (3), 241–246. Li, J., Chen, A., Yan, Z., Xu, G., Zhang, X., 1993. Particle–fluid contacting in circulating fluidized beds. AIChE, New York, pp. 49–54 (Preprint Volume for CFB-IV; Avidan, AA (Ed.)). Li, J., Ge, W., Wang, W., Yang, N., Liu, X., Wang, L., He, X., Wang, X., Wang, J., Kwauk, M., 2013. From Multiscale Modeling to Meso-science. Springer, New York. Li, J., Kwauk, M., 1994. Particle-Fluid Two-Phase Flow: The Energy-Minimization Multi-Scale Method. Metallurgical Industry Press, Beijing. Li, J., Tung, Y., Kwauk, M., 1988. Method of energy minimization in multi-scale modeling of particle–fluid two-phase flow. Circul. Fluid. Bed Technol. II, 89–103. Li, Y., Chen, B., Wang, F., Wang, Y., 1982. Hydrodynamic correlations for fast fluidization. In: First China–Japan Symposium of Fluidization, Hangzhou, China, pp. 124–134. Lim, K., Zhu, J., Grace, J., 1995. Hydrodynamics of gas–solid fluidization. Int. J. Multiph. Flow 21, 141–193. Liu, C., Zhao, M., Wang, W., Li, J., 2015. 3d cfd simulation of a circulating fluidized bed with on-line adjustment of mechanical valve. Chem. Eng. Sci. 137, 646–655. Liu, D., 1987. Set up the equations for tow-phase flows by the method of kinetic theory. Chin. J. Theor. Appl. Mech. 19 (3), 213–221 (In Chinese).

55

Liu, D., 1993. Fluid Dynamics of Two-Phase Systems. High Education Press, Beijing, China (In Chinese). Liu, L., 2011. Kinetic theory of aggregation in granular flow. AIChE J. 57 (12), 3331–3343. Manyele, S., Pärssinen, J., Zhu, J.-X., 2002. Characterizing particle aggregates in a high-density and high-flux cfb riser. Chem. Eng. J. 88 (1), 151–161. Marble, F.E., 1963. Dynamics of a gas containing small solid particles. In: Combustion and Propulsion (5th AGARDograph Colloquium). Pergamon Press, pp. 175– 213, New York. Mehrabadi, M., Tenneti, S., Garg, R., Subramaniam, S., 2015. Pseudo-turbulent gasphase velocity fluctuations in homogeneous gas–solid flow: fixed particle assemblies and freely evolving suspensions. J. Fluid Mech. 770, 210–246. Milioli, C.C., Milioli, F.E., Holloway, W., Agrawal, K., Sundaresan, S., 2013. Filtered two-fluid models of fluidized gas–particle flows: new constitutive relations. AIChE J. 59 (9), 3265–3275. OBrien, T.J., Syamlal, M., 1993. Particle cluster effects in the numerical simulation of a circulating fluidized bed. Circul. Fluid. Bed Technol. IV, 367–372. Ozel, A., Fede, P., Simonin, O., 2013. Development of filtered euler-euler two-phase model for circulating fluidised bed: high resolution simulation, formulation and a priori analyses. Int. J. Multiph. Flow 55, 43–63. Pai, S.-I., 1977. Mixture of fluid and solid particles. In: Two-Phase Flows. Springer, pp. 116–167, Berlin. Parmentier, J.-F., Simonin, O., Delsart, O., 2012. A functional subgrid drift velocity model for filtered drag prediction in dense fluidized bed. AIChE J. 58 (4), 1084–1098. Pärssinen, J., Zhu, J.-X., 2001a. Axial and radial solids distribution in a long and high-flux cfb riser. AIChE J. 47 (10), 2197–2205. Pärssinen, J., Zhu, J.-X., 2001b. Particle velocity and flow development in a long and high-flux circulating fluidized bed riser. Chem. Eng. Sci. 56 (18), 5295–5303. Pope, S.B., 2000. Turbulent Flows. Cambridge University Press, New York. Reh, L., 1971. Fluidized bed processing. Chem. Eng. Progr. 67 (2), 58–63. Reh, L., 2003. Development potentials and research needs in circulating fluidized bed combustion. China Particuol. 1 (5), 185–200. Schneiderbauer, S., Pirker, S., 2014. Filtered and heterogeneity-based subgrid modifications for gas–solid drag and solid stresses in bubbling fluidized beds. AIChE J. 60 (3), 839–854. Sharma, A.K., Tuzla, K., Matsen, J., Chen, J.C., 2000. Parametric effects of particle size and gas velocity on cluster characteristics in fast fluidized beds. Powder Technol. 111 (1), 114–122. Soong, C., Tuzla, K., Chen, J., 1994. Identification of particle clusters in circulating fluidized bed. Circul. Fluid. Bed Technol. 4, 615–620. Van der Hoef, M., van Sint Annaland, M., Deen, N., Kuipers, J., 2008. Numerical simulation of dense gas–solid fluidized beds: a multiscale modeling strategy. Annu. Rev. Fluid Mech. 40, 47–70. Wang, J., Ge, W., Li, J., 2008. Eulerian simulation of heterogeneous gas–solid flows in cfb risers: emms-based sub-grid scale model with a revised cluster description. Chem. Eng. Sci. 63 (6), 1553–1571. Wang, J., Zhou, Q., Hong, K., Wang, W., Li, J., 2012. An emms-based multi-fluid model (efm) for heterogeneous gas–solid riser flows: part ii. an alternative formulation from dominant mechanisms. Chem. Eng. Sci. 75, 349–358. Wang, W., Li, J., 2007. Simulation of gas–solid two-phase flow by a multi-scale cfd approach: extension of the emms model to the sub-grid scale level. Chem. Eng. Sci. 62 (1), 208–231. Xiao, H., Qi, H., You, C., Xu, X., 2003. Theoretical model of drag between gas and solid phase. J. Chem. Ind. Eng. (China) 54 (3), 311–315. Yang, N., Wang, W., Ge, W., Li, J., 2003. Cfd simulation of concurrent-up gas–solid flow in circulating fluidized beds with structure-dependent drag coefficient. Chem. Eng. J. 96 (1), 71–80. Yerushalmi, J., Turner, D.H., Squires, A.M., 1976. The fast fluidized bed. Ind. Eng. Chem. Process Des. Dev. 15 (1), 47–53. Zeneli, M., Nikolopoulos, A., Nikolopoulos, N., Grammelis, P., Kakaras, E., 2015. Application of an advanced coupled emms-tfm model to a pilot scale cfb carbonator. Chem. Eng. Sci. 138, 482–498. Zhao, B., Zhou, Q., Wang, J., Li, J., 2015. Cfd study of exit effect of high-density cfb risers with emms-based two-fluid model. Chem. Eng. Sci. 134, 477–488. Zhou, Q., Wang, J., 2015. Cfd study of mixing and segregation in cfb risers: extension of emms drag model to binary gas–solid flow. Chem. Eng. Sci. 122, 637–651. Zhou, Q., Wang, J., Li, J., 2014. Three-dimensional simulation of dense suspension upflow regime in high-density cfb risers with emms-based two-fluid model. Chem. Eng. Sci. 107, 206–217. Zhu, L.-T., Xie, L., Xiao, J., Luo, Z.-H., 2016. Filtered model for the cold-model gas– solid flow in a large-scale mto fluidized bed reactor. Chem. Eng. Sci. 143, 369–383. Zou, B., Li, H., Xia, Y., Ma, X., 1994. Cluster structure in a circulating fluidized bed. Powder Technol. 78 (2), 173–178.