An Assessment of Drag Models in Eulerian–Eulerian CFD Simulation of Gas–Solid Flow Hydrodynamics in Circulating Fluidized Bed Riser

: Accurate prediction of the hydrodynamic profile is important for circulating fluidized bed (CFB) reactor design and scale-up. Multiphase computational fluid dynamics (CFD) simulation with interphase momentum exchange is key to accurately predict the gas-solid profile along the height of the riser. The present work deals with the assessment of six different drag model capability to accurately predict the riser section axial solid holdup distribution in bench scale circulating fluidized bed. The difference between six drag model predictions were validated against the experiment data. Two-dimensional geometry, transient solver and Eulerian–Eulerian multiphase models were used. Six drag model simulation predictions were discussed with respect to axial and radial profile. The comparison between CFD simulation and experimental data shows that the Syamlal-O’Brien, Gidaspow, Wen-Yu and Huilin-Gidaspow drag models were successfully able to predict the riser upper section solid holdup distribution with better accuracy, however unable to predict the solid holdup transition region. On the other hand, the Gibilaro model and Helland drag model were successfully able to predict the bottom dense region, but the upper section solid holdup distribution was overpredicted. The CFD simulation comparison of different drag model has clearly shown the limitation of the drag model to accurately predict overall axial heterogeneity with accuracy. calculated from different height of the riser and averaged. CFD simulation results with the Helland drag model shows highest RMSE value compared to that Syamlal-O’Brien drag model has the lowest value. This comparison demonstrates that the Syamlal-O’Brien, Gidaspow and Wen-Yu drag models predict closer to the experimental data. Similarly, we found Gibilaro drag models RMSE value is higher than other drag model. Comparison between experiment and simulated axial solid holdup value using RMSE data demonstrates that the simulation prediction all along the height of riser reflected on RMSE value. Among all test drag models, none of the drag model successfully able to fully capture the upper dilute and bottom dense region. Therefore, in next part of this study is to formulate drag reduction correlation for given gas velocity, solid circulation rate and material type using three-dimensional (3D) CFD simulation.


Introduction
The gas-solid fluidized bed reactors are widely used in the petrochemical, power generation, environmental and metallurgical industries [1]. In the circulating fluidized bed (CFB) system, solid particles are separated from the fluid stream using cyclone and recycle to the bed. The solid particles fluid-like state provide excellent heat and mass transfer characteristic and circulation of solid offer operational flexibility [2,3]. The CFB riser section serves as the main reaction zone, in which strands of particles (clusters) influence the flow and the performance of reactor [4]. Over the years, considerable effort has been made to better design and operation of CFB reactor. In particular, focus to better understand the underlying gas-solid flow hydrodynamics in the riser section.
Over the year, computational fluid dynamics (CFD) domain has undergone significant development and regularly applied for investigating gas-solid multiphase flow phenomena. For example, direct numerical simulation (DNS), large eddy simulation (LES) and Reynolds-Averaged Navier-Stokes simulation (RANS) are three modelling methods at different computational scales to account the turbulence in multiphase flow. The accuracy of these numerical methods depends on the scheme. Several high order schemes were developed, such as essentially non-oscillatory (ENO) method [5,6], weighted ENO (WENO) method [7], discontinuous Galerkin (DG) method [8], radial basis function method [9] and gas-kinetic method [10,11]. However, for most industrial scale problem DNS and LES simulation were considered impractical for the general-purpose design tool [12]. RANS based CFD method has become a standard platform to simulate gas-solid multiphase flow in industrial scale fluidized bed reactor [13]. To model the gas-solid flow, the most commonly employed CFD simulation methods are the Eulerian-Eulerian two-fluid model (EE-TFM), [14,15,16] Eulerian-Lagrangian discrete element method (EL-DEM) [17,18], and hybrid multiphase particle-incell (MP-PIC) [19,20]. The principal difference among these approaches is on the treatment of the secondary solid phase. In EE-TFM approach kinetic theory of granular flow (KTGF) closure model were incorporated to simulate solid particle flow. EE-TFM model constitutive relationship contains solid phase stress, interphase momentum transfer and particle interactions [13,16,21]. Researchers have shown that the interphase drag force is the critical closure model to accurately predict the gassolid flows [22,23]. Agrawal et al. [23] has shown that the solid phase stress contribution was insignificant, and it is crucial to accurately account for the solid particle clusters.
Researchers have proposed the drag model based on flow conditions. In the literature, drag models were classically arranged into two groups, i.e., structure-based drag models and conventional drag models. Among the conventional drag models, Wen-Yu [24], Gidaspow [25] and Syamlal-O'Brien [26] were widely incorporated, which were derived from dense packed bed experimental pressure drop data and terminal velocity deduce using single solid particle. The most widely used Gidaspow model correlation consist of Ergun [27] and Wen-Yu [24] equations. To accurately simulate the riser gas-solid multiphase flow, researchers incorporated drag coefficient to account for the particle cluster. Researchers have pointed out that the particles cluster settling velocity is about 20% to 100% higher than the single particles. Also recognized that the such drag reduction occurred in the low solid volume fraction region. Various drag models have been incorporated to simulate the gassolid riser flow, as detail were provided in Appendix A. Among the works in the literature, Gidaspow and EMMS drag model has been widely used. Shah et al. [28] simulate fluid catalytic cracking (FCC) riser section using EMMS and the Gidaspow drag model for low and high solid circulation rate. They reported EMMS drag model successfully able to predict the riser axial heterogenous profile. However, unable to predict the solid volume fraction transition region. Vaishali et al. [29] perform CFD simulation for fast fluidization and dilute phase transport flow regime with Geldart B particles. They reported that the Wen-Yu drag model underpredicts the velocity profile and Syamlal-O'Brien drag model show better prediction. Wang et al. [30] investigate the models and model parameters variation on hydrodynamic prediction behavior for high solid circulation rate. Their drag model comparison study shows that the Syamlal-O'Brien, Gidaspow and Wen-Yu drag models predict the global profile similarly; however, on radial direction, the Syamlal-O'Brien drag model shows better accuracy. Similarly, Zhang et al. [31,32] compare EMMS and Gidaspow drag models, and it was reported that the EMMS drag model prediction is good in the axial and radial direction. Despite several studies, it is important to examine the influence of drag models for specific flow condition and material type. Hence, in the present work an assessment of the appropriate drag model selection for in-house CFB riser setup was carried out. Generated CFD simulation profile were plotted against experimentally obtain axial solid holdup data for detail assessment and validation purpose.
In the present work six drag models (Syamlal-O'Brien [26], Gidaspow [25], Wen-Yu [24], Huilin-Gidaspow [33], Gibilaro [34] and Helland [35]) were investigated. In order to keep the computational load at a reasonable level, simple 2D simulations were performed. A detailed discussion is presented to illustrate the inherent difference between different drag model prediction. Such a comparison study is expected to provide clear guidelines on selection of appropriate drag model selection for gas-solid riser flow.

Simulation System
In the present work, air used as fluidizing gas and 4 kg of silica sand solid particles were loaded into CFB reactor. CFD simulation was conducted at the same flow conditions as the work described in previous paper [16]. CFB riser simplified two-dimensional computational domain with the height of 3 m and a diameter of 0.025 m, shown in Figure 1a. To ensure symmetry flow the gas velocity was specified at the bottom, the solid phase were introduced from two sides and straight outlet, shown in Figure 1b. The gas and solid particle properties, system geometry and flow condition were summarized in Table 1. Our previous study investigated the optimum number of computational cells [16]. Therefore, 19,500 (13 × 1500) cells with a size of a 15-particle diameter were used to simulate the gas-solid flow in the riser.

CFD Model
Gas-solid CFD simulation was developed based on Eulerian-Eulerian two-fluid model (EE-TFM) [16] in ANSYS Fluent 2020R1 software. EE-TFM simulation equations are presented in Appendix B, and a detailed description of governing equations can be found elsewhere [14,15]. To correctly capture the riser gas-solid multiphase flow phenomena, granular temperature equation in partial-differential form was used to resemble the solid phase properties [36]. The standard k-ε turbulence model was selected to account the gas phase turbulence. Under Johnson and Jackson [16], the particle-wall interaction parameter, specularity coefficient (φ) and particle-wall restitution coefficient (ew) value were set at 0.0001 and 0.9, respectively. Elasticity between solid phase, defined as the particle-particle restitution coefficient (ess), was set to 0.9. The transient CFD calculations were performed, where time step of 0.0005 s was selected to capture riser gas-solid flow behavior. Coupling between the velocity and pressure was employed by a phase-coupled SIMPLE scheme, and other important simulation parameters are given in Table 2. Before post-processing, the riser outlet mass flow rate was monitored as a function of time to confirm the pseudo-steady state profile, as shown in Figure 2. All the transient simulations were run for 30 s and reported results were employed by time-averaging the simulation results over last 20 s. The interphase momentum exchange between the gas and solid phases was provided through the drag coefficient. Six different drag models investigated in the present work; drag correlations are given in Table 3.

Result and Discussion
Before comparing drag model simulation results against experimental data, riser outlet mass flow rate transient profile as function of time for the six different drag models were shown in Figure  3. Riser outlet mass flow rate data for all investigated drag models were plotted for about 20 s in which all simulations have achieved a pseudo-steady state. It is interesting to note that the outlet mass flow rate profile for the Syamlal-O'Brien [26], Gidaspow [25], Wen-Yu [24] and Huilin-Gidaspow [33] drag models were fluctuating around the average mass flow rate of 1.04 kg/s compared to the Gibilaro [34] and Helland [35] drag models. Such comparison shows a distinct mass flow rate profile for different drag models; it can also be helpful to interpret model prediction behavior by analyzing the outlet mass flow rate.    The influence of drag model is further illustrated by the time-averaged gas and solid phase axial velocity profile at a riser height of 1.5 m. The predicted time-averaged mean velocity for six different drag models were shown in Figure 5. It can be seen that the riser characteristic core-annulus flow structure observed in all cases. The highest gas velocity observed in the center and lowest value near the wall. Comparison clearly shows that the gas and the solid velocity magnitude and flow pattern variation for different drag models. For quantitative comparison, mean velocity value along the radial position were averaged and indicated in corresponding drag models. According to the value, gas phase axial average velocity is lowered by 4% from Syamlal-O'Brien, Gidaspow, Wen-Yu and Huilin-Gidaspow drag model to Gibilaro model and Helland drag model, whereas solid phase axial average velocity is lowered 39%. This over prediction of drag coefficient value results in higher mean solid velocity profile.
where experimental solid holdup value and solid holdup predicted from six different drag models. Small RMSE indicates a closer model prediction to the experimental data. Figure 7 displays the RMSE for CFD simulation prediction with respect to six different drag models. RMSE values were calculated from different height of the riser and averaged. CFD simulation results with the Helland drag model shows highest RMSE value compared to that Syamlal-O'Brien drag model has the lowest value. This comparison demonstrates that the Syamlal-O'Brien, Gidaspow and Wen-Yu drag models predict closer to the experimental data. Similarly, we found Gibilaro drag models RMSE value is higher than other drag model. Comparison between experiment and simulated axial solid holdup value using RMSE data demonstrates that the simulation prediction all along the height of riser reflected on RMSE value. Among all test drag models, none of the drag model successfully able to fully capture the upper dilute and bottom dense region. Therefore, in next part of this study is to formulate drag reduction correlation for given gas velocity, solid circulation rate and material type using three-dimensional (3D) CFD simulation.

Conclusions
In this work, we investigate the limitation of six different drag models to predict the solid holdup profile in the riser section of circulating fluidized bed (CFB) reactor. A classical two-fluid model (TFM) was used to simulate 2D computational domain in ANSYS Fluent 2020R1 software. The six different simulation results were compared against solid holdup data along the height of CFB riser. Comparison analysis shows that axial solid holdup predicted by the Syamlal-O'Brien [26], Gidaspow [25], Wen-Yu [24] and Huilin-Gidaspow [33] models are the closest to the experimental data above the height of 0.5 m. However, they underpredict the bottom solid holdup compared to the Helland and Gibilaro drag model. Major prediction differences were observed seen in the bottom section of the riser with Syamlal-O'Brien, Gidaspow, Wen-Yu and Huilin-Gidaspow drag models. A discrepancy between the six different simulation prediction was quantify using the root-mean-square error (RMSE) calculation; the Syamlal-O'Brien [26] drag model was found to be the most accurate, followed by the Wen-Yu and Gidaspow drag model. Overall, their performance shows the requirement of the incorporating modified drag coefficient for individual flow condition and material type to predict upper dilute and bottom dense region.  Table A1. Gas-solid drag models used in literature.

B.2 Momentum Conservation Equation
• Momentum equation: Solid and gas phase stress tensors: