- Email: [email protected]
Determination of choking in the EMMS model
Chemical Engineering Journal 357 (2019) 508–517
Contents lists available at ScienceDirect
Chemical Engineering Journal journal homepage: www.elsevier.com/locate/cej
Determination of choking in the EMMS model Jianhua Chen
T
State Key Laboratory of Multiphase Complex Systems, Institute of Process Engineering, Chinese Academy of Sciences, Beijing 100190, China
HIGHLIGHTS
GRAPHICAL ABSTRACT
the relationship between the • Elucidate two criteria of choking in the EMMS model.
choking criteria correspond to • Both the breakdown of the heterogeneous structures.
the GPU solver for the EMMS • Improve model to cope with the two criteria. a preferable computing con• Propose dition for the constrained FD regime. the impact of structural • Consolidate changes on choking via solution space visualization.
ARTICLE INFO
ABSTRACT
Keywords: Choking Fluidization EMMS model Regime transition Mesoscale
This paper investigates the relationship between the two choking criteria in the energy minimization multi-scale (EMMS) model. An improved GPU solver for the EMMS model is used to facilitate the related calculation. Both choking criteria agree well with the experiments in literature, and the predictions are almost identical to each other. The study confirms that both criteria define the critical conditions in terms of the breakdown of the fully heterogeneous structures. This structural change is due to the enhanced inter-phase interaction with increasing gas velocities at given solid flow rate. According to the numerical calculation, the solution space and the clustering constraint are visualized to consolidate the impact of structural changes on choking.
1. Introduction “Choking” is a characteristic phenomenon in the gas–solid fluidization and pneumatic conveying. The tremendous efforts on the choking research indicated its importance for understanding fluidization [1–4]. The works of Xu et al. [5,6] substantiated choking as the collapse of dilute suspension according to a thorough compilation of experimental data for saturation carrying capacity. Du et al. [7] investigated choking by examining real-time flow structures with the electrical capacitance tomography (ECT) technique. Their experiments revealed the structure variation around the transport velocity. Though various correlations have been proposed according to the massive experimental data, theoretical researches on choking are still scarce. The energy minimization multi-scale (EMMS) model [8] predicted choking as a sudden change in
voidage. The corresponding solid flow rate was commensurate with the so-called saturation carrying capacity for circulating fluidized beds (CFBs). In experiments, the S-shape axial solid concentration profile was typically observed in the riser [9], which indicated the coexistence of a dense phase at the bottom region and a dilute phase at the top region of the riser. Li et al. [10] pointed out that yielding the S-shape axial voidage was more reliable to identify choking than detecting the sudden change of flow structures. Then it was proven that a multi-scale computational fluid dynamics (CFD) approach coupling with the EMMS model [11,12] can quantitatively reproduce the bell-shape choking data [13]. Moreover, choking also has a significant impact on the phase transfer, because transfer coefficients usually depend on local heterogeneities such as clusters or bubbles to calculate non-dimensional numbers [14]. For instance, a so-called EMMS/mass algorithm [15,16] was proposed for CFD-
E-mail address: [email protected]. https://doi.org/10.1016/j.cej.2018.09.171 Received 11 June 2018; Received in revised form 22 August 2018; Accepted 21 September 2018 Available online 25 September 2018 1385-8947/ © 2018 Elsevier B.V. All rights reserved.
Chemical Engineering Journal 357 (2019) 508–517
J. Chen
Nomenclature
Ar CD dcl dp f g Gs Nst
Uc Uf Ug Upc Upf Upt Up Usc Usf Usi Wst
f max
Archimedes number, Ar = g ( p g ) gdp3 / µg2 , dimensionless drag coefficient, dimensionless cluster diameter, m particle diameter, m volume fraction of the dense phase, dimensionless gravity acceleration, m s 2 solid flow rate, kg m 2 s 1 energy consumption rate for suspending and transporting per unit mass of particles, W kg 1 superficial gas velocity in the dense phase, m s 1 superficial gas velocity in the dilute phase, m s 1 superficial gas velocity, m s 1 superficial solid velocity in the dense phase, m s 1 superficial solid velocity in the dilute phase, m s 1 superficial gas velocity at choking, m s 1 superficial solid velocity, m s 1 dense phase superficial slip velocity, m s 1 dilute phase superficial slip velocity, m s 1 inter-phase superficial slip velocity, m s 1 energy consumption rate with respect to unit volume of bed, W m 3
mf o uni
Abbreviations CFD Ck1 Ck2 CPU CUDA EMMS FD FDI GPU GRG PFC
c
computational fluid dynamics choking criterion-1 choking criterion-2 central processing unit compute unified device architecture energy minimization multi-scale fluid-dominated fluid-dominated ideal graphics processing unit generalized reduced gradient particle-fluid compromising
Subscripts cl c f g i p
Greek letters
µg
dilute phase voidage, dimensionless maximum voidage, dimensionless voidage at minimum fluidization, dimensionless c at local minimum of Nst , dimensionless voidage at uniform suspension, dimensionless
gas phase dynamic viscosity, Pa s density, kg m 3 average voidage, dimensionless dense phase voidage, dimensionless
cluster phase dense phase dilute phase gas phase inter-phase particle phase
coupled mass transfer computation. They found that the step change of structural parameters at choking can affect the simulating results remarkably. The variation of Sherwood number reported for CFB risers in literature was quantitatively explained owing to the choking phenomenon and clustering effects. In practice, two choking criteria have been proposed during the development of the EMMS model. The first one proposed by Li et al. [17] considers the preferred modes for the energy consumption rate at different regimes, where a particle–fluid compromising (PFC) regime and a fluid-dominated (FD) regime with heterogeneous structure are involved. This choking criterion is expressed as
(Wst )PFC = (Wst ) FD,
(1)
where Wst refers to the energy consumption rate for suspending and transporting particles with respect to unit volume of bed. Later on, Ge and Li [18] found the existence of a bi-stable characteristic with varying Nst versus c , that is, the second criterion
Nst (
mf )
= Nst ( o ),
(2)
Fig. 1. Comparison between model solvers and the experimental data taken from Li and Kwauk [13]. (FCC-air system with physical properties: particle density p = 929.5 kg/m3 ; gas density g = 1.1795 kg/m3 ; particle diameter d p = 54 µm ; gas dynamic viscosity µg = 1.8872 × 10 5 Pa·s ; average voidage at minimum fluidization mf = 0.5; solid flow rate Gs = 25 kg/(m2 s)).
where Nst refers to the energy consumption rate for suspending and transporting per unit mass of particles, o refers to a much higher value of the dense phase voidage of c which is more close to f rather than mf . It should be noted that o is obtained by the optimizing procedure instead of a constant value. When choking happens, there are two coexistent states with the same Nst value, but the related structural parameters are different. This criterion reveals that choking is an intrinsic hydrodynamic feature in view of the EMMS model. However, compared with the previous criterion based on regime analysis, the underlying physics remains vague. Since the choking prediction using these two criteria has been carried out by employing different numerical methods and software packages, no meaningful comparison on the performance of these criteria is available yet. In order to provide such a comparison, a unified GPU (Graphics Processing Unit) solver for the EMMS model is
developed and applied to choking prediction under various operating conditions. The rest of the paper is therefore arranged as follows: Section 2 introduces the numerical method to facilitate the analysis. Subsequently, Section 3 discusses the choking criterion-1 with a thorough study on the FD regime, while Section 4 reports the numerical characteristics of the choking criterion-2. Section 5 compares the two criteria according to the landscape of the solution space, which throws a light on further analytical works for choking. Finally the conclusions are drawn in Section 6. 509
Chemical Engineering Journal 357 (2019) 508–517
J. Chen
Fig. 2. The critical point of Ck1: (Wst )PFC = (Wst )FD , intersection of Wst curves pertaining to the regimes of PFC and FD respectively.
Fig. 4. Comparison of
2. Numerical method
uni
calculated by different methods.
method was developed and a rigorous numerical solution was obtained [18]. Since the traversal procedure costs heavy computing load, the GPU accelerated solver has been realized on the NVIDIA® CUDA platform with corresponding hardware, as illustrated in previous works [20,21].
The EMMS model is difficult to be solved analytically. In the budding period, the generalized reduced gradient (GRG) method and the nonlinear optimization software of GRG-2 were adopted to solve the model. Afterwards, in order to circumvent the possible missing of the global optimal solutions, Li et al. [19] developed a group of simplified submodels according to the GRG-2 numerical results. Then a direct traversal
2.1. GPU solver for the EMMS model There
are
eight
Fig. 3. Comparison of the structural parameters between PFC and FD. 510
independent
variables
X
=
Chemical Engineering Journal 357 (2019) 508–517
J. Chen
the GPU (NVIDIA® Tesla K80) program in comparison to more than 4 h for the CPU (Intel® E5-2680, 2.4 GHz) program. The solving speed is accelerated by several hundred times, which greatly facilitates the model application. For the three-dimensional discretization, computing time is prolonged by about 1.8 times due to more threads involved. Furthermore, this traversal algorithm accompanies GPU acceleration by an incident benefit of avoiding possible divergence, which can be extended to other nonlinear optimizing problems. The experimental data of FCC-air flow system [13] are used to validate the GPU solver, as shown in Fig. 1. The abrupt change in the voidage is quantitatively reproduced by the solver. Since both the experimental data and the EMMS model have been widely published, more detailed information on the experiments can refer to related references [13,19]. 3. Analysis on choking criterion-1 3.1. Structural parameters of PFC and FD Using the developed GPU solver for the EMMS model, the choking criterion-1 (Ck1, expressed by Eq. (1)) is investigated in the FCC-air system at a given solid flow rate of Gs = 25 kg/(m2 s), as shown in Fig. 2. The critical point of Upt = 2.28 m/s predicted by Ck1 locates at the intersection of the Wst curves pertaining to the regimes of PFC and FD respectively, where PFC refers to the stability condition of Nst = min with an additional constraint of c = mf , and FD refers to the stability condition of Nst = max with additional constraints of c = mf and f 0. It is noticeable that the extrema of Nst = max and Wst = min are equivalent to define FD according to the extremum characteristics [13]. Choosing Wst to compare PFC and FD is based on the assumption of least resistance in FD, so that the regime with smaller Wst will realize around choking. At lower Ug before this critical point, PFC prevails due to the lower Wst . While beyond choking, FD takes over PFC for the lower Wst . Furthermore, Li et al. [19] pointed out that beyond choking, another possible fluid-dominated ideal (FDI) regime of homogeneous suspension exists with an even lower Wst . The actual regime for Ug > Upt would lie somewhere between FD and FDI. Fortunately, this uncertainty of the actual flow state is not influential to the intersection characteristic of the Wst curves. However, the exact meaning of the FD regime in Ck1 has not been well understood yet. As a first step, the structural parameters are compared between PFC and FD, as shown in Fig. 3. Given the constraint of c = mf in both regimes, the discrepancies of structural parameters such as f and are not large, whereas a qualitatively difference exists in f , as shown in Fig. 3(a)-(b). It should be mentioned that the closed values of f and between PFC and FD do not mean insignificant effects of the stability condition on fluidization status, because at given operating conditions, only one regime of PFC and FD can realize while the other is a pseudo state. Moreover, the additional constraints exert considerable influence on the condition of Nst = max. If the unconstrained Nst = max is used as the stability condition, a homogeneous suspension state will come out, as shown in Fig. 3(a). Comparison between Fig. 3(c) and (d) indicates a major difference in velocities, that is, PFC corresponds to a near-zero Upf while FD corresponds to near-zero values of both Upc and Uc . In other words, the solid flow rate is mainly contributed from the dense phase in PFC but from the dilute phase in FD. Obviously, the two regimes demand different energy consumption rates except at the critical condition. In addition, the curves for Uc and Upc nearly overlap for the close values in both regimes, indicating a much higher drag coefficient within the dense phase than the globally averaged drag coefficient, which is supported by the study of Li and Kwauk [22]. The strange behaviors of structural parameters call for careful deliberation of the FD regime.
Fig. 5. Bi-stable characteristic of choking criterion-2.
[ f , c, f , Uf , Uc, Upf , Upc, d cl ] in the EMMS model. The five conservation equations combined with the cluster diameter correlation provide six constraints. Detailed formulas and definitions are shown in Appendix A. From the viewpoint of traversal solving, two variables can be selected as the solution space coordinates. They are discretized into sufficiently small steps according to the physical meaning. Then the objective function of Nst is calculated in the solution space. Since the calculation in each traversal node is completely independent, the computing is fairly suitable to be parallelized on the GPU. That is, taking [ f , c ] as the solution space coordinate, both f and c are divided into 3072 intervals. To perform the discretization on the GPU kernel, the threads are conducted into thread blocks and thread grids hierarchically on the CUDA platform as Block dimension: dim3 dimBlock (16, 16); Grid dimension: dim3 dimGrid (3072, 12); CPU sends these coordinates to GPU and the conservation equations are solved on each GPU thread. All Nst values satisfying the constraints are recorded and collected. Then the minimum Nst is searched directly on the CPU once the solution space is obtained. The detailed solving scheme is shown in Appendix B. CUDA extends the programing language by allowing the definition of C functions denoted by GPU kernels. When the kernel is called, the functions are executed by individual threads. Inasmuch as the threads are wrapped into blocks and grids, the traversal algorithm is implemented simultaneously. Though looping on [ f , c ] seems enough for solving the EMMS model, an internal iteration on f is indispensable due to the non-independence between the inter-phase force balance equation and the cluster correlation, i.e. Eqs. (A.5) and (A.6) respectively, which leads to a complicated solution pattern in the two-dimensional space [18]. Alternatively, it is feasible to traverse over three variables by including f, thereby forming a three-dimensional solution space. In comparison to the two-dimensional discretization in the previous works [20,21], three-dimensional resolution of [ f , c, f ] is also performed to explore the landscape of the solution space. In this case, the discretization of 256 × 256 × 256 is utilized and the threads are conducted as Block dimension: dim3 dimBlock (16, 16); Grid dimension: dim3 dimGrid (256, 256). 2.2. Validation of the solver The GPU acceleration effects can be explained by the two-dimensional discretization case of FCC-air system at the operating condition of Ug = 2.0 m/s, Gs = 25 kg/(m2 s). The computing time costs 42 s for 511
Chemical Engineering Journal 357 (2019) 508–517
J. Chen
Fig. 6. Numerical characteristics related to choking criterion-2.
3.2. The equivalence between
f
=
uni
and Uc
0 in FD
velocity us at uniform state equals to the terminal velocity of single particle ut ,
0 Hopefully, the obsession of FD along with f = uni and Uc [23,19] may be surmounted owing to the improved solver. In fact, there 0 in FD. In order to explain is an equivalence between f = uni and Uc this equivalence, the definition of uni should be assessed at first. There are at least two ways to calculate uni . First one is directly from the relationship between the solid flux and the gas flux, that is uni,1
=1
Up Ug
.
us = Ug /
uni
Up/(1
uni )
ut ,
(4)
and hence uni,2
= [(Ug + Up + u t )
(Ug + Up + ut )2 4ut Ug ]/2ut.
(5)
The calculated uni and f in FD are shown in Fig. 4, all of which are in good agreement with each other, especially at higher gas velocities. In the current solver, indeed, the stability condition of Nst = max| c= mf for FD always comes out a homogeneous state with f = 0. The additional
(3)
Another one is according to Matsen [24] that the two-phase slip 512
Chemical Engineering Journal 357 (2019) 508–517
J. Chen
Ug = fUc + (1 f ) Uf
(1 f ) Uf ,
Up = fUpc + (1 f ) Upf
(1 f ) Upf .
(6) (7)
Then it is reasonable to use Eq. (3) to calculate the dilute phase voidage, regarding Eqs. 6,7, there is f
=1
Upf Uf
1
Up Ug
=
uni.
(8)
0 is almost equivalent to f = uni , noticing that uni Therefore, Uc itself is not necessarily precise. On the other hand, if f = uni is speci0 will certainly come out for the FD regime. It fied, the solution of Uc 0 should be mentioned that the equivalence between f = uni and Uc may be invalid if the constraint of c = mf is not satisfied when the cluster structure completely breaks down, because Uc and Upc would not approach to zero simultaneously thereafter. 4. Analysis on choking criterion-2 Fig. 7. Comparison of predicted saturation carrying capacities with the experimental data in references [13,26]. (Physical properties: FCC system is the same as that in Fig. 1; Chalmers system: particle density p = 2600 kg/m3 ; gas density g = 0.32 kg/m3 ; particle diameter d p = 260 µm ; gas dynamic viscosity µg = 4.832× 10 5 Pa·s ; average voidage at minimum fluidization mf = 0.43).
Fig. 8. Cluster constraints (contours of (Gs = 25 kg/(m2 s)).
The choking criterion-2 (Ck2, expressed by Eq. (2)) resorts to a rigorous numerical calculation with the unconstrained stability condition of Nst = min . With increasing Ug , there may have two local minima for Nst , one of which corresponds to c = mf , and the other corresponds to c = o , as shown in Fig. 5. At the critical condition, the bi-stable characteristic with identical Nst indicates the regime transition as a shift between different structural parameters. Fig. 6 shows the numerical characteristics of the structural parameters and the energy consumption rates. As shown in Fig. 6(a), at the critical point which corresponds to the jump change in , c suddenly increases from mf to o , and f tends to 0.5 to guarantee the momentum balance after jump change. However, there is no sudden change in f , which means the abrupt change in the dense phase has not influenced the dilute phase remarkably. Fig. 6(b) shows the velocity profiles for the dense phase and the dilute phase. Uc and Upc separate with each other and Uc approaches to Uf beyond the jump change. The solid flux is still mainly contributed by the dense phase beyond choking, concerning the difference between Upc and Upf . If the solutions of both Nst ( mf ) and Nst ( o ) are recorded, the corresponding curves of Nst versus Ug intersect at the critical point, as shown in Fig. 6(c). In other words, the curve of Nst vs. Ug driven by Nst = min is continuous, but its derivative is discontinuous [25]. As for Wst calculated by Nst = min, it exhibits an abrupt change around Ug = 2.28 m/s, declining to an even lower value than that for the homogeneous suspension, Wst = gGs , as shown in Fig. 6(d). However, since the relationship Wst = Nst × p (1 ) holds, the sharp decrease in Wst is compensated by the jump change in , so that a continuous Nst retains at the two critical points. Noticed that the cluster diameter of dcl is a reciprocal function of Nst , it is therefore continuous at choking, as seen in Fig. 6(e). In addition, Fig. 6(f) shows the superficial slip velocities. The calculated Usf tends to be the terminal velocity of ut which is about 0.08 m/s in this case, and Usc is much smaller than Usf , even with the sharp increase beyond choking. Most interestingly, the averaged slip velocity of Us seems equal to Usi before choking, which means the slip velocity is mainly determined by the inter-phase interaction in this regime. At choking, both Us and Usi decline to be small values comparable to Usf , indicating that the gas–solid flow changes from strong-slip to weak-slip modes. At present, it remains unsettled of the actual state for regimes after the jump change. Fortunately, in the sense of engineering application, both FDI and Nst = min can reproduce similar average voidage with an acceptable accuracy in this regime.
Usi = 0) in the solution space
0 as well as an artificial lower bound for f is therefore constraint of f indispensable to reproduce FD. Nonethelss, if Uc = 0 is specified along with c = mf , the solver reproduces the FD curve opportunely, even though the condition of Nst = max is replaced by Wst = min . In this sense, Wst = min|[ c= mf , Uc= 0] is also a feasible computing condition for FD, and it is more easily to be implemented in comparison to the foregoing constraint. Indeed, the condition of [ c = mf , f = uni] neglecting extremum of Nst is enough to reproduce FD in the analytical works [23,19], the reason of which will be explained in Section 5.1. 0 according to the weak slip 0 leads to Upc Consequently, Uc velocity in the dense phase constrained by c = mf . At this condition, we have
5. Comparison between the two choking criteria Fig. 7 shows the predicted saturation carrying capacities at choking by Ck1 and Ck2. The calculated results are compared with the experimental data in references [13,26], and detailed physical properties can refer to Ge and Li [18]. The comparison confirms that both criteria are 513
Chemical Engineering Journal 357 (2019) 508–517
J. Chen
Fig. 9. Relationship between the contours of Usi = 0 and Wst at different Ug : (a) Ug = 1.0 m/s; (b) Ug = 2.28 m/s.
operating conditions. As aforementioned, the FD solution locates on the slice of f = uni which is rendered by the contour of Nst . It is clearly seen that the intersection among c = mf , f = uni and Usi = 0 determines the FD solution. It further expounds why the analytical work [23,19] can reproduce the FD curve with given c and f neglecting the extremum of Nst . 5.2. Further discussion These researches verified that two choking criteria reach the same goal by different means, i.e., they identify choking by intersection of Wst and Nst curves, respectively. Essentially, choking is the critical condition that the aggregate structures cannot retain any longer due to the intense interaction between the dilute phase and the dense phase. The change of either f for Ck1 or c for Ck2 indicates another more preferred mode, leading to an inevitable regime transition. The first choking criterion Ck1 was reacquainted with a thorough 0 and understanding of the FD regime. The equivalence between Uc f = uni plays an important role to determine the structural parameters. 0 in FD means that the solid flux is almost conFurthermore, Uc tributed from the dilute phase in comparison to the dense phase in PFC, which indicates an alternative preponderance between the dilute phase and the dense phase in FD and PFC, respectively. The physical idea of Ck1 has already been discussed in the light of mesoscience [27]. That is, the transition from A-B compromising regime to B-dominated regime can be defined by a critical condition of Bmin |A B = BA B , where Bmin |A B refers to the quantity of B under the condition of B = min however constrained by a characterizing structure of the regime of A-B compromising, whilst BA B refers to the quantity of B given by the regime of A-B compromising. In the choking case, correspondingly, the (Wst )min |[ c = mf , Uc= 0] = (Wst ) Nst = min , criterion becomes taking [ c = mf , Uc = 0] as the characterizing structure of the regime of A-B compromising, i.e. Nst = min . In other words, imaging the resident 0 ) exerted by a bypassing homogeneous dense clusters ( c mf , Uc suspension ( f uni ), when the energy consumption is lower than that 1), the transition of the distinct heterogeneous state ( c mf , f occurs because the input power is capable to break up the clusters. Present work exemplifies this profound inference by clarifying the condition for calculating FD. It should be mentioned that this hypothesis for defining regime transition was inspired by physical intuition rather than rigorous mathematical proofs. As a preliminary strategy, it deserves further work to be more reliable when extended to other systems.
Fig. 10. Calculated choking lines with varying parameters.
feasible to predict choking. It is obvious that only slight difference exists between the predicted choking lines by Ck1 and Ck2, which indicates an inner link between the two choking criteria. Indeed, the clustering correlation imposes a strong constraint upon the solution space, so that the mesoscale structure is crucial to determine the real solution. Therefore, a visualization of the solution space is helpful to understand this impact. 5.1. Landscape of the solution space According to the previous study [25], the EMMS model solution should satisfy the constraint exerted by the cluster diameter equation, which means all feasible solutions including the choking point must locate on the iso-surface of Usi = 0 in the solution space, where Usi refers to the superficial slip velocity between the dilute phase and the dense phase, and the definition of Usi can refer to Appendix B. It should be mentioned that Usi is non-zero in spite of Usi = 0, as illustrated in Fig. 6(f). Fig. 8 shows three typical contours of Usi = 0 at the operating condition of Gs = 25 kg/(m2 s). It is obvious the contour of Usi = 0 evolves in the solution space with varying gas velocities. Moreover, Fig. 9 illustrates the relationship between the contours of Usi = 0 and Wst for the solution of PFC and FD under different 514
Chemical Engineering Journal 357 (2019) 508–517
J. Chen
The second choking criterion Ck2 was proposed according to the bistable characteristics of Nst . Although it benefits choking prediction with an unconstrained stability condition, it should be cautious to apply Nst = min to the FD regime due to several uncertainties. For instance, the dilute phase voidage of f does not decrease after the abrupt change, though particles are prone to penetrate into the dilute phase due to the collapse of clusters. The cluster correlation may be a main limitation, because it is not necessarily valid once the dense phase structure has been destructed beyond choking. As pointed out by Li et al. [28], it is still unsolved for theoretical cluster modeling due to its dynamic characteristics, which leaves a weakness of the EMMS model. Despite the existing uncertainties, current model and software provide a feasible tool to predict the flow structures and choking for gas-solid fluidization. Moreover, it is capable to reflect the influence of physical properties on choking. This feature enables a possible extension to bi-disperse fluidization systems [29,30] which are commonly encountered in industrial applications. In the bi-disperse fluidization, different particle inertia and particle–particle interactions account for the modeling challenges. To this end, complicated kinetic theory models [31,32] and modified drag models [33,34] have been incorporated into the CFD platforms to cope with the binary mixture. Hopefully, current framework can provide a preliminary estimation of global flow parameters, in particular, choking for bi-disperse systems. Taking the aforementioned FCC-air system as an example, if we increase the Archimedes number Ar by nearly 10 times via adjusting one physical parameter and keeping others unchanged, the calculated choking lines will move with reasonable trends, that is, the choking velocity augments with increasing d p or decreasing µg , and it declines with increasing g , as seen in Fig. 10. These results qualitatively coincide with the experimental research on poly-disperse fluidization [35]. Their work found that the superficial gas velocity at choking will increase with higher percentage of larger particles. Nonetheless, it is noteworthy that the EMMS model is developed for mono-disperse
systems. Therefore, a detailed extension to bi-disperse fluidization requires further work and the particle-particle interaction should be taken into account. 6. Conclusions The structural characteristics of the two criteria in the EMMS model were investigated respectively owing to the improved GPU solver. The study on the first criterion revealed an equivalence between f = uni 0 in the FD regime, according to which a modified computing and Uc condition for FD was proposed. We confirmed the two choking criteria in the EMMS model to be essentially identical, in terms of the breakdown of the fully heterogeneous structures with increasing Ug or decreasing Gs due to the enhanced inter-phase interaction. In general, both criteria in the EMMS model were proven to work well compared with the experimental data in literature. The dependence on clustering constraints indicates that mesoscale structures are essential to understand choking in fluidization. Further works on modeling mesoscale structures are required to circumvent these system-specific and dynamic complexities. Acknowledgements This work is financially supported by the State Key Laboratory of Multiphase Complex Systems (Grant No. MPCS-2015-A-03) and the Key Research Program of Frontier Science (Grant No. QYZDJ-SSW-JSC029), as well as the Strategic Priority Research Program (Grant No. XDA21030700), Chinese Academy of Sciences. The author is grateful to Prof. Jinghai Li for his instruction and encouragement of this research, and Prof. Wei Ge for his pertinent and constructive suggestions. The colleagues of the EMMS group and their helpful discussions are also greatly appreciated.
Appendix A. Formulas of the EMMS model See Table 1.
Table 1 Conservation equations, the stability condition and correlations. Equation meanings Continuity of the gas phase
Continuity of the particle phase
Cluster diameter correlation
dcl =
Stability condition
(A.2)
=
1 1
c
=(
f 1
Nst = Ug
c Upc 1 c
Rec =
Drag coefficient of single particle
g dp Usc µg
CDc0 =
24 Rec
+
Drag coefficient of particle swarms
CDc = CDc0
c
(A.4)
g) g
(A.5)
3.6 Rec0.313 4.7
g dp Usf µg
CDf0 =
24 Re f
CDf = CDf0 Nst,mf = dp Umf g µg
515
2
+
+ f
g) g p
g) g p
p (1
Ref =
( p
( p
(A.6) (A.7)
= min
(A.8)
)
(
f Upf 1 f
Usf = Uf
1.75 3 mf
(A.3)
g) g
p
p
f (1 f ) Uf
Wst = Nst ×
Usc = Uc
(
2 g Usi 3 C = ( c)( p g ) g 4 Di d cl ( p g ) gUp / [ p (1 max )] Nst,mf dp Nst Nst,mf
Volumetric energy consumption rate
Equation for calculating Umf
(A.1)
2 g Usf 3 CDf 4 dp
Force balance of the inter-phase
Equation for calculating Nst,mf
Ug = fUc + (1 f ) Uf 2 g Usc 3 CDc 4 dp
Force balance of the dilute phase
Reynolds number
Eqs.
Up = fUpc + (1 f ) Upf
Force balance of the dense phase
Superficial slip velocity
Formulas
Usi = Uf
3.6 Re 0.313 f 4.7
Umf +
f Upc 1 c
Rei =
g d cl Usi µg
CDi0 =
24 Rei
+
CDi = CDi0 (1 Up mf 1 mf
150(1 mf ) dp Umf g 3 µg mf
=
dp3 g ( p µg2
g) g
)
(1 f )
3.6 Rei0.313 f ) 4.7
(A.9) (A.10) (A.11) (A.12) (A.13) (A.14)
Chemical Engineering Journal 357 (2019) 508–517
J. Chen
Appendix B. The traversal solving scheme See Fig. A1
Fig. A1.
flows, Chem. Eng. Sci. 116 (2014) 710–724. [15] W. Dong, W. Wang, J. Li, A multiscale mass transfer model for gas-solid riser flows: Part I—sub-grid model and simple tests, Chem. Eng. Sci. 63 (2008) 2798–2810. [16] W. Dong, W. Wang, J. Li, A multiscale mass transfer model for gas-solid riser flows: Part II—sub-grid simulation of ozone decomposition, Chem. Eng. Sci. 63 (2008) 2811–2823. [17] J. Li, M. Kwauk, L. Reh, Role of energy minimization in gas/solid fluidization, in: O. E. Potter, D. Nicklin (Eds.), Fluidization 7th, Engineering Foundation, 1992, p. 83. [18] W. Ge, J. Li, Physical mapping of fluidization regimes—the EMMS approach, Chem. Eng. Sci. 57 (2002) 3993–4004. [19] J. Li, C. Cheng, Z. Zhang, J. Yuan, A. Nemet, F.N. Fett, The EMMS model — its application, development and updated concepts, Chem. Eng. Sci. 54 (1999) 5409–5425. [20] W. Ge, W. Wang, N. Yang, J. Li, M. Kwauk, F. Chen, J. Chen, X. Fang, L. Guo, X. He, X. Liu, Y. Liu, B. Lu, J. Wang, J. Wang, L. Wang, X. Wang, Q. Xiong, M. Xu, L. Deng, Y. Han, C. Hou, L. Hua, W. Huang, B. Li, C. Li, F. Li, Y. Ren, J. Xu, N. Zhang, Y. Zhang, G. Zhou, G. Zhou, Meso-scale oriented simulation towards virtual process engineering (VPE) — the EMMS paradigm, Chem. Eng. Sci. 66 (2011) 4426–4458. [21] M. Du, S. Hu, J. Chen, X. Liu, W. Ge, Extremum characteristics of energy consumption in fluidization analyzed by using EMMS, Chem. Eng. J. 342 (2018) 386–394. [22] J. Li, M. Kwauk, Multiscale nature of complex fluid-particle systems, Ind. Eng. Chem. Res. 40 (2001) 4227–4237. [23] G. Xu, J. Li, Analytical solution of the energy-minimization multi-scale model for gas-solid two-phase flow, Chem. Eng. Sci. 53 (1998) 1349–1366. [24] J. Matsen, Mechanism of choking and entrainment, Powder Technol. 32 (1982) 21–23. [25] J.H. Chen, N. Yang, W. Ge, J.H. Li, Stability-driven structure evolution: exploring the intrinsic similarity between gas-solid and gas-liquid systems, Chin. J. Chem. Eng. 20 (2012) 167–177. [26] M. Couturier, B. Doucette, D. Stevens, Temperature, gas concentration and solid
References [1] H.T. Bi, J.R. Grace, J.X. Zhu, Types of choking in vertical pneumatic systems, Int. J. Multiphase Flow 19 (1993) 1077–1092. [2] H. Bi, J. Zhu, Static instability analysis of circulating fluidized beds and concept of high-density risers, AIChE J. 39 (1993) 1272–1280. [3] D. Bai, A.S. Issangya, J.R. Grace, A novel method for determination of choking velocities, Powder Technol. 97 (1998) 59–62. [4] W.-C. Yang, Choking revisited, Ind. Eng. Chem. Res. 43 (2004) 5496–5506. [5] G. Xu, K. Nomura, S. Gao, K. Kato, More fundamentals of dilute suspension collapse and choking for vertical conveying systems, AIChE J. 47 (2001) 2177–2196. [6] G. Xu, S. Gao, Dependence of the solid flow rate on gas velocity at choking of different types, Powder Technol. 150 (2005) 185–188. [7] B. Du, W. Warsito, L.-S. Fan, ECT studies of the choking phenomenon in a gas-solid circulating fluidized bed, AIChE J. 50 (2004) 1386–1406. [8] J. Li, Y. Tung, M. Kwauk, Method of energy minimization in multi-scale modeling of particle-fluid two-phase flow, in: L.J. Basu (Ed.), Circulating Fluidized Bed Technology II, Pergamon Press, New York, 1988, pp. 75–89. [9] Y. Li, M. Kwauk, The Dynamics of Fast Fluidization, Springer, US, Boston, MA, 1980, pp. 537–544. [10] J. Li, L. Wen, W. Ge, H. Cui, J. Ren, Dissipative structure in concurrent-up gas-solid flow, Chem. Eng. Sci. 53 (1998) 3367–3379. [11] W. Wang, B. Lu, J. Li, Choking and flow regime transitions: Simulation by a multiscale CFD approach, Chem. Eng. Sci. 62 (2007) 814–819. [12] W. Wang, B. Lu, N. Zhang, Z. Shi, J. Li, A review of multiscale CFD for gas-solid CFB modeling, Int. J. Multiph. Flow 36 (2010) 109–118. [13] J. Li, M. Kwauk, Particle-Fluid Two-Phase Flow—The Energy-Minimization MultiScale Method, Metallurgical Industry Press, Beijing, 1994. [14] N.G. Deen, E.A.J.F. Peters, J.T. Padding, J.A.M. Kuipers, Review of direct numerical simulation of fluid-particle mass, momentum and heat transfer in dense gas-solid
516
Chemical Engineering Journal 357 (2019) 508–517
J. Chen
[27] [28] [29] [30]
mass flux profiles within a large circulating fluidized bed combustor, in: J.F. Large (Ed.), Proceedings of the 12th international conference on fluidized bed combustion, ASME, New York, 1991, pp. 107–116. J. Li, W. Huang, From multiscale to mesoscience: Addressing mesoscales in mesoregimes of different levels, Annu. Rev. Chem. Biomolecular Eng. 9 (2018) 41–60. J. Li, W. Ge, W. Wang, N. Yang, W. Huang, Focusing on mesoscales: from the energy-minimization multiscale model to mesoscience, Curr. Opin. Chem. Eng. 13 (2016) 10–23. O. Owoyemi, L. Mazzei, P. Lettieri, CFD modeling of binary-fluidized suspensions and investigation of role of particle-particle drag on mixing and segregation, AIChE J. 53 (2007) 1924–1940. J.D. Wilde, G. Richards, S. Benyahia, Qualitative numerical study of simultaneous high-G-intensified gas-solids contact, separation and segregation in a bi-disperse
[31] [32] [33] [34] [35]
517
rotating fluidized bed in a vortex chamber, Adv. Powder Technol. 27 (2016) 1453–1463. H. Iddir, H. Arastoopour, Modeling of multitype particle flow using the kinetic theory approach, AIChE J. 51 (2005) 1620–1632. S. Benyahia, Verification and validation study of some polydisperse kinetic theories, Chem. Eng. Sci. 63 (2008) 5672–5680. Q. Zhou, J. Wang, CFD study of mixing and segregation in CFB risers: extension of EMMS drag model to binary gas-solid flow, Chem. Eng. Sci. 122 (2015) 637–651. M. Lungu, Y. Zhou, J. Wang, Y. Yang, A CFD study of a bi-disperse gas-solid fluidized bed: Effect of the EMMS sub grid drag correction, Powder Technol. 280 (2015) 154–172. H. Arastoopour, S.C. Lin, S.A. Weil, Analysis of vertical pneumatic conveying of solids using multiphase flow models, AIChE J. 28 (1982) 467–473.
Contents lists available at ScienceDirect
Chemical Engineering Journal journal homepage: www.elsevier.com/locate/cej
Determination of choking in the EMMS model Jianhua Chen
T
State Key Laboratory of Multiphase Complex Systems, Institute of Process Engineering, Chinese Academy of Sciences, Beijing 100190, China
HIGHLIGHTS
GRAPHICAL ABSTRACT
the relationship between the • Elucidate two criteria of choking in the EMMS model.
choking criteria correspond to • Both the breakdown of the heterogeneous structures.
the GPU solver for the EMMS • Improve model to cope with the two criteria. a preferable computing con• Propose dition for the constrained FD regime. the impact of structural • Consolidate changes on choking via solution space visualization.
ARTICLE INFO
ABSTRACT
Keywords: Choking Fluidization EMMS model Regime transition Mesoscale
This paper investigates the relationship between the two choking criteria in the energy minimization multi-scale (EMMS) model. An improved GPU solver for the EMMS model is used to facilitate the related calculation. Both choking criteria agree well with the experiments in literature, and the predictions are almost identical to each other. The study confirms that both criteria define the critical conditions in terms of the breakdown of the fully heterogeneous structures. This structural change is due to the enhanced inter-phase interaction with increasing gas velocities at given solid flow rate. According to the numerical calculation, the solution space and the clustering constraint are visualized to consolidate the impact of structural changes on choking.
1. Introduction “Choking” is a characteristic phenomenon in the gas–solid fluidization and pneumatic conveying. The tremendous efforts on the choking research indicated its importance for understanding fluidization [1–4]. The works of Xu et al. [5,6] substantiated choking as the collapse of dilute suspension according to a thorough compilation of experimental data for saturation carrying capacity. Du et al. [7] investigated choking by examining real-time flow structures with the electrical capacitance tomography (ECT) technique. Their experiments revealed the structure variation around the transport velocity. Though various correlations have been proposed according to the massive experimental data, theoretical researches on choking are still scarce. The energy minimization multi-scale (EMMS) model [8] predicted choking as a sudden change in
voidage. The corresponding solid flow rate was commensurate with the so-called saturation carrying capacity for circulating fluidized beds (CFBs). In experiments, the S-shape axial solid concentration profile was typically observed in the riser [9], which indicated the coexistence of a dense phase at the bottom region and a dilute phase at the top region of the riser. Li et al. [10] pointed out that yielding the S-shape axial voidage was more reliable to identify choking than detecting the sudden change of flow structures. Then it was proven that a multi-scale computational fluid dynamics (CFD) approach coupling with the EMMS model [11,12] can quantitatively reproduce the bell-shape choking data [13]. Moreover, choking also has a significant impact on the phase transfer, because transfer coefficients usually depend on local heterogeneities such as clusters or bubbles to calculate non-dimensional numbers [14]. For instance, a so-called EMMS/mass algorithm [15,16] was proposed for CFD-
E-mail address: [email protected]. https://doi.org/10.1016/j.cej.2018.09.171 Received 11 June 2018; Received in revised form 22 August 2018; Accepted 21 September 2018 Available online 25 September 2018 1385-8947/ © 2018 Elsevier B.V. All rights reserved.
Chemical Engineering Journal 357 (2019) 508–517
J. Chen
Nomenclature
Ar CD dcl dp f g Gs Nst
Uc Uf Ug Upc Upf Upt Up Usc Usf Usi Wst
f max
Archimedes number, Ar = g ( p g ) gdp3 / µg2 , dimensionless drag coefficient, dimensionless cluster diameter, m particle diameter, m volume fraction of the dense phase, dimensionless gravity acceleration, m s 2 solid flow rate, kg m 2 s 1 energy consumption rate for suspending and transporting per unit mass of particles, W kg 1 superficial gas velocity in the dense phase, m s 1 superficial gas velocity in the dilute phase, m s 1 superficial gas velocity, m s 1 superficial solid velocity in the dense phase, m s 1 superficial solid velocity in the dilute phase, m s 1 superficial gas velocity at choking, m s 1 superficial solid velocity, m s 1 dense phase superficial slip velocity, m s 1 dilute phase superficial slip velocity, m s 1 inter-phase superficial slip velocity, m s 1 energy consumption rate with respect to unit volume of bed, W m 3
mf o uni
Abbreviations CFD Ck1 Ck2 CPU CUDA EMMS FD FDI GPU GRG PFC
c
computational fluid dynamics choking criterion-1 choking criterion-2 central processing unit compute unified device architecture energy minimization multi-scale fluid-dominated fluid-dominated ideal graphics processing unit generalized reduced gradient particle-fluid compromising
Subscripts cl c f g i p
Greek letters
µg
dilute phase voidage, dimensionless maximum voidage, dimensionless voidage at minimum fluidization, dimensionless c at local minimum of Nst , dimensionless voidage at uniform suspension, dimensionless
gas phase dynamic viscosity, Pa s density, kg m 3 average voidage, dimensionless dense phase voidage, dimensionless
cluster phase dense phase dilute phase gas phase inter-phase particle phase
coupled mass transfer computation. They found that the step change of structural parameters at choking can affect the simulating results remarkably. The variation of Sherwood number reported for CFB risers in literature was quantitatively explained owing to the choking phenomenon and clustering effects. In practice, two choking criteria have been proposed during the development of the EMMS model. The first one proposed by Li et al. [17] considers the preferred modes for the energy consumption rate at different regimes, where a particle–fluid compromising (PFC) regime and a fluid-dominated (FD) regime with heterogeneous structure are involved. This choking criterion is expressed as
(Wst )PFC = (Wst ) FD,
(1)
where Wst refers to the energy consumption rate for suspending and transporting particles with respect to unit volume of bed. Later on, Ge and Li [18] found the existence of a bi-stable characteristic with varying Nst versus c , that is, the second criterion
Nst (
mf )
= Nst ( o ),
(2)
Fig. 1. Comparison between model solvers and the experimental data taken from Li and Kwauk [13]. (FCC-air system with physical properties: particle density p = 929.5 kg/m3 ; gas density g = 1.1795 kg/m3 ; particle diameter d p = 54 µm ; gas dynamic viscosity µg = 1.8872 × 10 5 Pa·s ; average voidage at minimum fluidization mf = 0.5; solid flow rate Gs = 25 kg/(m2 s)).
where Nst refers to the energy consumption rate for suspending and transporting per unit mass of particles, o refers to a much higher value of the dense phase voidage of c which is more close to f rather than mf . It should be noted that o is obtained by the optimizing procedure instead of a constant value. When choking happens, there are two coexistent states with the same Nst value, but the related structural parameters are different. This criterion reveals that choking is an intrinsic hydrodynamic feature in view of the EMMS model. However, compared with the previous criterion based on regime analysis, the underlying physics remains vague. Since the choking prediction using these two criteria has been carried out by employing different numerical methods and software packages, no meaningful comparison on the performance of these criteria is available yet. In order to provide such a comparison, a unified GPU (Graphics Processing Unit) solver for the EMMS model is
developed and applied to choking prediction under various operating conditions. The rest of the paper is therefore arranged as follows: Section 2 introduces the numerical method to facilitate the analysis. Subsequently, Section 3 discusses the choking criterion-1 with a thorough study on the FD regime, while Section 4 reports the numerical characteristics of the choking criterion-2. Section 5 compares the two criteria according to the landscape of the solution space, which throws a light on further analytical works for choking. Finally the conclusions are drawn in Section 6. 509
Chemical Engineering Journal 357 (2019) 508–517
J. Chen
Fig. 2. The critical point of Ck1: (Wst )PFC = (Wst )FD , intersection of Wst curves pertaining to the regimes of PFC and FD respectively.
Fig. 4. Comparison of
2. Numerical method
uni
calculated by different methods.
method was developed and a rigorous numerical solution was obtained [18]. Since the traversal procedure costs heavy computing load, the GPU accelerated solver has been realized on the NVIDIA® CUDA platform with corresponding hardware, as illustrated in previous works [20,21].
The EMMS model is difficult to be solved analytically. In the budding period, the generalized reduced gradient (GRG) method and the nonlinear optimization software of GRG-2 were adopted to solve the model. Afterwards, in order to circumvent the possible missing of the global optimal solutions, Li et al. [19] developed a group of simplified submodels according to the GRG-2 numerical results. Then a direct traversal
2.1. GPU solver for the EMMS model There
are
eight
Fig. 3. Comparison of the structural parameters between PFC and FD. 510
independent
variables
X
=
Chemical Engineering Journal 357 (2019) 508–517
J. Chen
the GPU (NVIDIA® Tesla K80) program in comparison to more than 4 h for the CPU (Intel® E5-2680, 2.4 GHz) program. The solving speed is accelerated by several hundred times, which greatly facilitates the model application. For the three-dimensional discretization, computing time is prolonged by about 1.8 times due to more threads involved. Furthermore, this traversal algorithm accompanies GPU acceleration by an incident benefit of avoiding possible divergence, which can be extended to other nonlinear optimizing problems. The experimental data of FCC-air flow system [13] are used to validate the GPU solver, as shown in Fig. 1. The abrupt change in the voidage is quantitatively reproduced by the solver. Since both the experimental data and the EMMS model have been widely published, more detailed information on the experiments can refer to related references [13,19]. 3. Analysis on choking criterion-1 3.1. Structural parameters of PFC and FD Using the developed GPU solver for the EMMS model, the choking criterion-1 (Ck1, expressed by Eq. (1)) is investigated in the FCC-air system at a given solid flow rate of Gs = 25 kg/(m2 s), as shown in Fig. 2. The critical point of Upt = 2.28 m/s predicted by Ck1 locates at the intersection of the Wst curves pertaining to the regimes of PFC and FD respectively, where PFC refers to the stability condition of Nst = min with an additional constraint of c = mf , and FD refers to the stability condition of Nst = max with additional constraints of c = mf and f 0. It is noticeable that the extrema of Nst = max and Wst = min are equivalent to define FD according to the extremum characteristics [13]. Choosing Wst to compare PFC and FD is based on the assumption of least resistance in FD, so that the regime with smaller Wst will realize around choking. At lower Ug before this critical point, PFC prevails due to the lower Wst . While beyond choking, FD takes over PFC for the lower Wst . Furthermore, Li et al. [19] pointed out that beyond choking, another possible fluid-dominated ideal (FDI) regime of homogeneous suspension exists with an even lower Wst . The actual regime for Ug > Upt would lie somewhere between FD and FDI. Fortunately, this uncertainty of the actual flow state is not influential to the intersection characteristic of the Wst curves. However, the exact meaning of the FD regime in Ck1 has not been well understood yet. As a first step, the structural parameters are compared between PFC and FD, as shown in Fig. 3. Given the constraint of c = mf in both regimes, the discrepancies of structural parameters such as f and are not large, whereas a qualitatively difference exists in f , as shown in Fig. 3(a)-(b). It should be mentioned that the closed values of f and between PFC and FD do not mean insignificant effects of the stability condition on fluidization status, because at given operating conditions, only one regime of PFC and FD can realize while the other is a pseudo state. Moreover, the additional constraints exert considerable influence on the condition of Nst = max. If the unconstrained Nst = max is used as the stability condition, a homogeneous suspension state will come out, as shown in Fig. 3(a). Comparison between Fig. 3(c) and (d) indicates a major difference in velocities, that is, PFC corresponds to a near-zero Upf while FD corresponds to near-zero values of both Upc and Uc . In other words, the solid flow rate is mainly contributed from the dense phase in PFC but from the dilute phase in FD. Obviously, the two regimes demand different energy consumption rates except at the critical condition. In addition, the curves for Uc and Upc nearly overlap for the close values in both regimes, indicating a much higher drag coefficient within the dense phase than the globally averaged drag coefficient, which is supported by the study of Li and Kwauk [22]. The strange behaviors of structural parameters call for careful deliberation of the FD regime.
Fig. 5. Bi-stable characteristic of choking criterion-2.
[ f , c, f , Uf , Uc, Upf , Upc, d cl ] in the EMMS model. The five conservation equations combined with the cluster diameter correlation provide six constraints. Detailed formulas and definitions are shown in Appendix A. From the viewpoint of traversal solving, two variables can be selected as the solution space coordinates. They are discretized into sufficiently small steps according to the physical meaning. Then the objective function of Nst is calculated in the solution space. Since the calculation in each traversal node is completely independent, the computing is fairly suitable to be parallelized on the GPU. That is, taking [ f , c ] as the solution space coordinate, both f and c are divided into 3072 intervals. To perform the discretization on the GPU kernel, the threads are conducted into thread blocks and thread grids hierarchically on the CUDA platform as Block dimension: dim3 dimBlock (16, 16); Grid dimension: dim3 dimGrid (3072, 12); CPU sends these coordinates to GPU and the conservation equations are solved on each GPU thread. All Nst values satisfying the constraints are recorded and collected. Then the minimum Nst is searched directly on the CPU once the solution space is obtained. The detailed solving scheme is shown in Appendix B. CUDA extends the programing language by allowing the definition of C functions denoted by GPU kernels. When the kernel is called, the functions are executed by individual threads. Inasmuch as the threads are wrapped into blocks and grids, the traversal algorithm is implemented simultaneously. Though looping on [ f , c ] seems enough for solving the EMMS model, an internal iteration on f is indispensable due to the non-independence between the inter-phase force balance equation and the cluster correlation, i.e. Eqs. (A.5) and (A.6) respectively, which leads to a complicated solution pattern in the two-dimensional space [18]. Alternatively, it is feasible to traverse over three variables by including f, thereby forming a three-dimensional solution space. In comparison to the two-dimensional discretization in the previous works [20,21], three-dimensional resolution of [ f , c, f ] is also performed to explore the landscape of the solution space. In this case, the discretization of 256 × 256 × 256 is utilized and the threads are conducted as Block dimension: dim3 dimBlock (16, 16); Grid dimension: dim3 dimGrid (256, 256). 2.2. Validation of the solver The GPU acceleration effects can be explained by the two-dimensional discretization case of FCC-air system at the operating condition of Ug = 2.0 m/s, Gs = 25 kg/(m2 s). The computing time costs 42 s for 511
Chemical Engineering Journal 357 (2019) 508–517
J. Chen
Fig. 6. Numerical characteristics related to choking criterion-2.
3.2. The equivalence between
f
=
uni
and Uc
0 in FD
velocity us at uniform state equals to the terminal velocity of single particle ut ,
0 Hopefully, the obsession of FD along with f = uni and Uc [23,19] may be surmounted owing to the improved solver. In fact, there 0 in FD. In order to explain is an equivalence between f = uni and Uc this equivalence, the definition of uni should be assessed at first. There are at least two ways to calculate uni . First one is directly from the relationship between the solid flux and the gas flux, that is uni,1
=1
Up Ug
.
us = Ug /
uni
Up/(1
uni )
ut ,
(4)
and hence uni,2
= [(Ug + Up + u t )
(Ug + Up + ut )2 4ut Ug ]/2ut.
(5)
The calculated uni and f in FD are shown in Fig. 4, all of which are in good agreement with each other, especially at higher gas velocities. In the current solver, indeed, the stability condition of Nst = max| c= mf for FD always comes out a homogeneous state with f = 0. The additional
(3)
Another one is according to Matsen [24] that the two-phase slip 512
Chemical Engineering Journal 357 (2019) 508–517
J. Chen
Ug = fUc + (1 f ) Uf
(1 f ) Uf ,
Up = fUpc + (1 f ) Upf
(1 f ) Upf .
(6) (7)
Then it is reasonable to use Eq. (3) to calculate the dilute phase voidage, regarding Eqs. 6,7, there is f
=1
Upf Uf
1
Up Ug
=
uni.
(8)
0 is almost equivalent to f = uni , noticing that uni Therefore, Uc itself is not necessarily precise. On the other hand, if f = uni is speci0 will certainly come out for the FD regime. It fied, the solution of Uc 0 should be mentioned that the equivalence between f = uni and Uc may be invalid if the constraint of c = mf is not satisfied when the cluster structure completely breaks down, because Uc and Upc would not approach to zero simultaneously thereafter. 4. Analysis on choking criterion-2 Fig. 7. Comparison of predicted saturation carrying capacities with the experimental data in references [13,26]. (Physical properties: FCC system is the same as that in Fig. 1; Chalmers system: particle density p = 2600 kg/m3 ; gas density g = 0.32 kg/m3 ; particle diameter d p = 260 µm ; gas dynamic viscosity µg = 4.832× 10 5 Pa·s ; average voidage at minimum fluidization mf = 0.43).
Fig. 8. Cluster constraints (contours of (Gs = 25 kg/(m2 s)).
The choking criterion-2 (Ck2, expressed by Eq. (2)) resorts to a rigorous numerical calculation with the unconstrained stability condition of Nst = min . With increasing Ug , there may have two local minima for Nst , one of which corresponds to c = mf , and the other corresponds to c = o , as shown in Fig. 5. At the critical condition, the bi-stable characteristic with identical Nst indicates the regime transition as a shift between different structural parameters. Fig. 6 shows the numerical characteristics of the structural parameters and the energy consumption rates. As shown in Fig. 6(a), at the critical point which corresponds to the jump change in , c suddenly increases from mf to o , and f tends to 0.5 to guarantee the momentum balance after jump change. However, there is no sudden change in f , which means the abrupt change in the dense phase has not influenced the dilute phase remarkably. Fig. 6(b) shows the velocity profiles for the dense phase and the dilute phase. Uc and Upc separate with each other and Uc approaches to Uf beyond the jump change. The solid flux is still mainly contributed by the dense phase beyond choking, concerning the difference between Upc and Upf . If the solutions of both Nst ( mf ) and Nst ( o ) are recorded, the corresponding curves of Nst versus Ug intersect at the critical point, as shown in Fig. 6(c). In other words, the curve of Nst vs. Ug driven by Nst = min is continuous, but its derivative is discontinuous [25]. As for Wst calculated by Nst = min, it exhibits an abrupt change around Ug = 2.28 m/s, declining to an even lower value than that for the homogeneous suspension, Wst = gGs , as shown in Fig. 6(d). However, since the relationship Wst = Nst × p (1 ) holds, the sharp decrease in Wst is compensated by the jump change in , so that a continuous Nst retains at the two critical points. Noticed that the cluster diameter of dcl is a reciprocal function of Nst , it is therefore continuous at choking, as seen in Fig. 6(e). In addition, Fig. 6(f) shows the superficial slip velocities. The calculated Usf tends to be the terminal velocity of ut which is about 0.08 m/s in this case, and Usc is much smaller than Usf , even with the sharp increase beyond choking. Most interestingly, the averaged slip velocity of Us seems equal to Usi before choking, which means the slip velocity is mainly determined by the inter-phase interaction in this regime. At choking, both Us and Usi decline to be small values comparable to Usf , indicating that the gas–solid flow changes from strong-slip to weak-slip modes. At present, it remains unsettled of the actual state for regimes after the jump change. Fortunately, in the sense of engineering application, both FDI and Nst = min can reproduce similar average voidage with an acceptable accuracy in this regime.
Usi = 0) in the solution space
0 as well as an artificial lower bound for f is therefore constraint of f indispensable to reproduce FD. Nonethelss, if Uc = 0 is specified along with c = mf , the solver reproduces the FD curve opportunely, even though the condition of Nst = max is replaced by Wst = min . In this sense, Wst = min|[ c= mf , Uc= 0] is also a feasible computing condition for FD, and it is more easily to be implemented in comparison to the foregoing constraint. Indeed, the condition of [ c = mf , f = uni] neglecting extremum of Nst is enough to reproduce FD in the analytical works [23,19], the reason of which will be explained in Section 5.1. 0 according to the weak slip 0 leads to Upc Consequently, Uc velocity in the dense phase constrained by c = mf . At this condition, we have
5. Comparison between the two choking criteria Fig. 7 shows the predicted saturation carrying capacities at choking by Ck1 and Ck2. The calculated results are compared with the experimental data in references [13,26], and detailed physical properties can refer to Ge and Li [18]. The comparison confirms that both criteria are 513
Chemical Engineering Journal 357 (2019) 508–517
J. Chen
Fig. 9. Relationship between the contours of Usi = 0 and Wst at different Ug : (a) Ug = 1.0 m/s; (b) Ug = 2.28 m/s.
operating conditions. As aforementioned, the FD solution locates on the slice of f = uni which is rendered by the contour of Nst . It is clearly seen that the intersection among c = mf , f = uni and Usi = 0 determines the FD solution. It further expounds why the analytical work [23,19] can reproduce the FD curve with given c and f neglecting the extremum of Nst . 5.2. Further discussion These researches verified that two choking criteria reach the same goal by different means, i.e., they identify choking by intersection of Wst and Nst curves, respectively. Essentially, choking is the critical condition that the aggregate structures cannot retain any longer due to the intense interaction between the dilute phase and the dense phase. The change of either f for Ck1 or c for Ck2 indicates another more preferred mode, leading to an inevitable regime transition. The first choking criterion Ck1 was reacquainted with a thorough 0 and understanding of the FD regime. The equivalence between Uc f = uni plays an important role to determine the structural parameters. 0 in FD means that the solid flux is almost conFurthermore, Uc tributed from the dilute phase in comparison to the dense phase in PFC, which indicates an alternative preponderance between the dilute phase and the dense phase in FD and PFC, respectively. The physical idea of Ck1 has already been discussed in the light of mesoscience [27]. That is, the transition from A-B compromising regime to B-dominated regime can be defined by a critical condition of Bmin |A B = BA B , where Bmin |A B refers to the quantity of B under the condition of B = min however constrained by a characterizing structure of the regime of A-B compromising, whilst BA B refers to the quantity of B given by the regime of A-B compromising. In the choking case, correspondingly, the (Wst )min |[ c = mf , Uc= 0] = (Wst ) Nst = min , criterion becomes taking [ c = mf , Uc = 0] as the characterizing structure of the regime of A-B compromising, i.e. Nst = min . In other words, imaging the resident 0 ) exerted by a bypassing homogeneous dense clusters ( c mf , Uc suspension ( f uni ), when the energy consumption is lower than that 1), the transition of the distinct heterogeneous state ( c mf , f occurs because the input power is capable to break up the clusters. Present work exemplifies this profound inference by clarifying the condition for calculating FD. It should be mentioned that this hypothesis for defining regime transition was inspired by physical intuition rather than rigorous mathematical proofs. As a preliminary strategy, it deserves further work to be more reliable when extended to other systems.
Fig. 10. Calculated choking lines with varying parameters.
feasible to predict choking. It is obvious that only slight difference exists between the predicted choking lines by Ck1 and Ck2, which indicates an inner link between the two choking criteria. Indeed, the clustering correlation imposes a strong constraint upon the solution space, so that the mesoscale structure is crucial to determine the real solution. Therefore, a visualization of the solution space is helpful to understand this impact. 5.1. Landscape of the solution space According to the previous study [25], the EMMS model solution should satisfy the constraint exerted by the cluster diameter equation, which means all feasible solutions including the choking point must locate on the iso-surface of Usi = 0 in the solution space, where Usi refers to the superficial slip velocity between the dilute phase and the dense phase, and the definition of Usi can refer to Appendix B. It should be mentioned that Usi is non-zero in spite of Usi = 0, as illustrated in Fig. 6(f). Fig. 8 shows three typical contours of Usi = 0 at the operating condition of Gs = 25 kg/(m2 s). It is obvious the contour of Usi = 0 evolves in the solution space with varying gas velocities. Moreover, Fig. 9 illustrates the relationship between the contours of Usi = 0 and Wst for the solution of PFC and FD under different 514
Chemical Engineering Journal 357 (2019) 508–517
J. Chen
The second choking criterion Ck2 was proposed according to the bistable characteristics of Nst . Although it benefits choking prediction with an unconstrained stability condition, it should be cautious to apply Nst = min to the FD regime due to several uncertainties. For instance, the dilute phase voidage of f does not decrease after the abrupt change, though particles are prone to penetrate into the dilute phase due to the collapse of clusters. The cluster correlation may be a main limitation, because it is not necessarily valid once the dense phase structure has been destructed beyond choking. As pointed out by Li et al. [28], it is still unsolved for theoretical cluster modeling due to its dynamic characteristics, which leaves a weakness of the EMMS model. Despite the existing uncertainties, current model and software provide a feasible tool to predict the flow structures and choking for gas-solid fluidization. Moreover, it is capable to reflect the influence of physical properties on choking. This feature enables a possible extension to bi-disperse fluidization systems [29,30] which are commonly encountered in industrial applications. In the bi-disperse fluidization, different particle inertia and particle–particle interactions account for the modeling challenges. To this end, complicated kinetic theory models [31,32] and modified drag models [33,34] have been incorporated into the CFD platforms to cope with the binary mixture. Hopefully, current framework can provide a preliminary estimation of global flow parameters, in particular, choking for bi-disperse systems. Taking the aforementioned FCC-air system as an example, if we increase the Archimedes number Ar by nearly 10 times via adjusting one physical parameter and keeping others unchanged, the calculated choking lines will move with reasonable trends, that is, the choking velocity augments with increasing d p or decreasing µg , and it declines with increasing g , as seen in Fig. 10. These results qualitatively coincide with the experimental research on poly-disperse fluidization [35]. Their work found that the superficial gas velocity at choking will increase with higher percentage of larger particles. Nonetheless, it is noteworthy that the EMMS model is developed for mono-disperse
systems. Therefore, a detailed extension to bi-disperse fluidization requires further work and the particle-particle interaction should be taken into account. 6. Conclusions The structural characteristics of the two criteria in the EMMS model were investigated respectively owing to the improved GPU solver. The study on the first criterion revealed an equivalence between f = uni 0 in the FD regime, according to which a modified computing and Uc condition for FD was proposed. We confirmed the two choking criteria in the EMMS model to be essentially identical, in terms of the breakdown of the fully heterogeneous structures with increasing Ug or decreasing Gs due to the enhanced inter-phase interaction. In general, both criteria in the EMMS model were proven to work well compared with the experimental data in literature. The dependence on clustering constraints indicates that mesoscale structures are essential to understand choking in fluidization. Further works on modeling mesoscale structures are required to circumvent these system-specific and dynamic complexities. Acknowledgements This work is financially supported by the State Key Laboratory of Multiphase Complex Systems (Grant No. MPCS-2015-A-03) and the Key Research Program of Frontier Science (Grant No. QYZDJ-SSW-JSC029), as well as the Strategic Priority Research Program (Grant No. XDA21030700), Chinese Academy of Sciences. The author is grateful to Prof. Jinghai Li for his instruction and encouragement of this research, and Prof. Wei Ge for his pertinent and constructive suggestions. The colleagues of the EMMS group and their helpful discussions are also greatly appreciated.
Appendix A. Formulas of the EMMS model See Table 1.
Table 1 Conservation equations, the stability condition and correlations. Equation meanings Continuity of the gas phase
Continuity of the particle phase
Cluster diameter correlation
dcl =
Stability condition
(A.2)
=
1 1
c
=(
f 1
Nst = Ug
c Upc 1 c
Rec =
Drag coefficient of single particle
g dp Usc µg
CDc0 =
24 Rec
+
Drag coefficient of particle swarms
CDc = CDc0
c
(A.4)
g) g
(A.5)
3.6 Rec0.313 4.7
g dp Usf µg
CDf0 =
24 Re f
CDf = CDf0 Nst,mf = dp Umf g µg
515
2
+
+ f
g) g p
g) g p
p (1
Ref =
( p
( p
(A.6) (A.7)
= min
(A.8)
)
(
f Upf 1 f
Usf = Uf
1.75 3 mf
(A.3)
g) g
p
p
f (1 f ) Uf
Wst = Nst ×
Usc = Uc
(
2 g Usi 3 C = ( c)( p g ) g 4 Di d cl ( p g ) gUp / [ p (1 max )] Nst,mf dp Nst Nst,mf
Volumetric energy consumption rate
Equation for calculating Umf
(A.1)
2 g Usf 3 CDf 4 dp
Force balance of the inter-phase
Equation for calculating Nst,mf
Ug = fUc + (1 f ) Uf 2 g Usc 3 CDc 4 dp
Force balance of the dilute phase
Reynolds number
Eqs.
Up = fUpc + (1 f ) Upf
Force balance of the dense phase
Superficial slip velocity
Formulas
Usi = Uf
3.6 Re 0.313 f 4.7
Umf +
f Upc 1 c
Rei =
g d cl Usi µg
CDi0 =
24 Rei
+
CDi = CDi0 (1 Up mf 1 mf
150(1 mf ) dp Umf g 3 µg mf
=
dp3 g ( p µg2
g) g
)
(1 f )
3.6 Rei0.313 f ) 4.7
(A.9) (A.10) (A.11) (A.12) (A.13) (A.14)
Chemical Engineering Journal 357 (2019) 508–517
J. Chen
Appendix B. The traversal solving scheme See Fig. A1
Fig. A1.
flows, Chem. Eng. Sci. 116 (2014) 710–724. [15] W. Dong, W. Wang, J. Li, A multiscale mass transfer model for gas-solid riser flows: Part I—sub-grid model and simple tests, Chem. Eng. Sci. 63 (2008) 2798–2810. [16] W. Dong, W. Wang, J. Li, A multiscale mass transfer model for gas-solid riser flows: Part II—sub-grid simulation of ozone decomposition, Chem. Eng. Sci. 63 (2008) 2811–2823. [17] J. Li, M. Kwauk, L. Reh, Role of energy minimization in gas/solid fluidization, in: O. E. Potter, D. Nicklin (Eds.), Fluidization 7th, Engineering Foundation, 1992, p. 83. [18] W. Ge, J. Li, Physical mapping of fluidization regimes—the EMMS approach, Chem. Eng. Sci. 57 (2002) 3993–4004. [19] J. Li, C. Cheng, Z. Zhang, J. Yuan, A. Nemet, F.N. Fett, The EMMS model — its application, development and updated concepts, Chem. Eng. Sci. 54 (1999) 5409–5425. [20] W. Ge, W. Wang, N. Yang, J. Li, M. Kwauk, F. Chen, J. Chen, X. Fang, L. Guo, X. He, X. Liu, Y. Liu, B. Lu, J. Wang, J. Wang, L. Wang, X. Wang, Q. Xiong, M. Xu, L. Deng, Y. Han, C. Hou, L. Hua, W. Huang, B. Li, C. Li, F. Li, Y. Ren, J. Xu, N. Zhang, Y. Zhang, G. Zhou, G. Zhou, Meso-scale oriented simulation towards virtual process engineering (VPE) — the EMMS paradigm, Chem. Eng. Sci. 66 (2011) 4426–4458. [21] M. Du, S. Hu, J. Chen, X. Liu, W. Ge, Extremum characteristics of energy consumption in fluidization analyzed by using EMMS, Chem. Eng. J. 342 (2018) 386–394. [22] J. Li, M. Kwauk, Multiscale nature of complex fluid-particle systems, Ind. Eng. Chem. Res. 40 (2001) 4227–4237. [23] G. Xu, J. Li, Analytical solution of the energy-minimization multi-scale model for gas-solid two-phase flow, Chem. Eng. Sci. 53 (1998) 1349–1366. [24] J. Matsen, Mechanism of choking and entrainment, Powder Technol. 32 (1982) 21–23. [25] J.H. Chen, N. Yang, W. Ge, J.H. Li, Stability-driven structure evolution: exploring the intrinsic similarity between gas-solid and gas-liquid systems, Chin. J. Chem. Eng. 20 (2012) 167–177. [26] M. Couturier, B. Doucette, D. Stevens, Temperature, gas concentration and solid
References [1] H.T. Bi, J.R. Grace, J.X. Zhu, Types of choking in vertical pneumatic systems, Int. J. Multiphase Flow 19 (1993) 1077–1092. [2] H. Bi, J. Zhu, Static instability analysis of circulating fluidized beds and concept of high-density risers, AIChE J. 39 (1993) 1272–1280. [3] D. Bai, A.S. Issangya, J.R. Grace, A novel method for determination of choking velocities, Powder Technol. 97 (1998) 59–62. [4] W.-C. Yang, Choking revisited, Ind. Eng. Chem. Res. 43 (2004) 5496–5506. [5] G. Xu, K. Nomura, S. Gao, K. Kato, More fundamentals of dilute suspension collapse and choking for vertical conveying systems, AIChE J. 47 (2001) 2177–2196. [6] G. Xu, S. Gao, Dependence of the solid flow rate on gas velocity at choking of different types, Powder Technol. 150 (2005) 185–188. [7] B. Du, W. Warsito, L.-S. Fan, ECT studies of the choking phenomenon in a gas-solid circulating fluidized bed, AIChE J. 50 (2004) 1386–1406. [8] J. Li, Y. Tung, M. Kwauk, Method of energy minimization in multi-scale modeling of particle-fluid two-phase flow, in: L.J. Basu (Ed.), Circulating Fluidized Bed Technology II, Pergamon Press, New York, 1988, pp. 75–89. [9] Y. Li, M. Kwauk, The Dynamics of Fast Fluidization, Springer, US, Boston, MA, 1980, pp. 537–544. [10] J. Li, L. Wen, W. Ge, H. Cui, J. Ren, Dissipative structure in concurrent-up gas-solid flow, Chem. Eng. Sci. 53 (1998) 3367–3379. [11] W. Wang, B. Lu, J. Li, Choking and flow regime transitions: Simulation by a multiscale CFD approach, Chem. Eng. Sci. 62 (2007) 814–819. [12] W. Wang, B. Lu, N. Zhang, Z. Shi, J. Li, A review of multiscale CFD for gas-solid CFB modeling, Int. J. Multiph. Flow 36 (2010) 109–118. [13] J. Li, M. Kwauk, Particle-Fluid Two-Phase Flow—The Energy-Minimization MultiScale Method, Metallurgical Industry Press, Beijing, 1994. [14] N.G. Deen, E.A.J.F. Peters, J.T. Padding, J.A.M. Kuipers, Review of direct numerical simulation of fluid-particle mass, momentum and heat transfer in dense gas-solid
516
Chemical Engineering Journal 357 (2019) 508–517
J. Chen
[27] [28] [29] [30]
mass flux profiles within a large circulating fluidized bed combustor, in: J.F. Large (Ed.), Proceedings of the 12th international conference on fluidized bed combustion, ASME, New York, 1991, pp. 107–116. J. Li, W. Huang, From multiscale to mesoscience: Addressing mesoscales in mesoregimes of different levels, Annu. Rev. Chem. Biomolecular Eng. 9 (2018) 41–60. J. Li, W. Ge, W. Wang, N. Yang, W. Huang, Focusing on mesoscales: from the energy-minimization multiscale model to mesoscience, Curr. Opin. Chem. Eng. 13 (2016) 10–23. O. Owoyemi, L. Mazzei, P. Lettieri, CFD modeling of binary-fluidized suspensions and investigation of role of particle-particle drag on mixing and segregation, AIChE J. 53 (2007) 1924–1940. J.D. Wilde, G. Richards, S. Benyahia, Qualitative numerical study of simultaneous high-G-intensified gas-solids contact, separation and segregation in a bi-disperse
[31] [32] [33] [34] [35]
517
rotating fluidized bed in a vortex chamber, Adv. Powder Technol. 27 (2016) 1453–1463. H. Iddir, H. Arastoopour, Modeling of multitype particle flow using the kinetic theory approach, AIChE J. 51 (2005) 1620–1632. S. Benyahia, Verification and validation study of some polydisperse kinetic theories, Chem. Eng. Sci. 63 (2008) 5672–5680. Q. Zhou, J. Wang, CFD study of mixing and segregation in CFB risers: extension of EMMS drag model to binary gas-solid flow, Chem. Eng. Sci. 122 (2015) 637–651. M. Lungu, Y. Zhou, J. Wang, Y. Yang, A CFD study of a bi-disperse gas-solid fluidized bed: Effect of the EMMS sub grid drag correction, Powder Technol. 280 (2015) 154–172. H. Arastoopour, S.C. Lin, S.A. Weil, Analysis of vertical pneumatic conveying of solids using multiphase flow models, AIChE J. 28 (1982) 467–473.