Fluidized beds modeling: Validation of 2D and 3D simulations against experiments

Accepted Manuscript Fluidized beds modeling: Validation of 2D and 3D simulations against experiments

Hongbo Shi, Alexandra Komrakova, Petr Nikrityuk PII: DOI: Reference:

S0032-5910(18)30947-1 https://doi.org/10.1016/j.powtec.2018.11.043 PTEC 13878

To appear in:

Powder Technology

Received date: Revised date: Accepted date:

13 December 2017 20 October 2018 10 November 2018

Please cite this article as: Hongbo Shi, Alexandra Komrakova, Petr Nikrityuk , Fluidized beds modeling: Validation of 2D and 3D simulations against experiments. Ptec (2018), https://doi.org/10.1016/j.powtec.2018.11.043

This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

ACCEPTED MANUSCRIPT Fluidized Beds Modeling: Validation of 2D and 3D Simulations Against Experiments Hongbo Shi1, Alexandra Komrakova1*, Petr Nikrityuk2 1

Department of Mechanical Engineering, University of Alberta, Edmonton, Alberta, Canada

T

2

IP

Department of Chemical and Materials Engineering, University of Alberta,

Abstract

CR

Edmonton, Alberta, Canada

US

The flow field in fluidized beds is studied numerically using a two-phase Eulerian-Eulerian model. The main goal of the study is to demonstrate the importance of the numerical model

AN

setup: we investigate the effect of model dimensionality (2D & 3D), the impact of the flow regimes (laminar & turbulent), the choice of the model parameters (specularity coefficient

M

responsible for particle-wall interaction and restitution coefficient characterizing particle-particle interaction), mesh resolution, and gas-solid drag sub-models

ED

(Syamlal-O’Brien and Gidaspow). The results of simulations are validated with the experimental data obtained by Taghipour et. al. (Chem. Eng. Sci., 60(24), 2005,

PT

6857-6867). The comparison of the results shows that the 2D simulations tend to overpredict the experimental data in terms of void fraction distribution across the flow.

CE

When the 3D simulations are used, the numerically obtained profiles of the void fraction inside the bed agree well with experimentally measured values. The results of the 3D

AC

simulations clearly illustrate that the finite size of the domain where fluidization occurs plays a significant role in the size of bubbles and in the spatial distribution of the volume fraction of particles inside fluidization zone. It is shown that turbulence has a minor influence on the fluidized bed structure and void fraction distribution. However, it significantly affects the gas phase velocity in the upstream region of a reactor. A thorough study of numerous numerical model setups concludes that for the 2D simulations of fluidized beds there is a combination of model parameters such as specularity and

*

Corresponding author. E-mail: [email protected]

ACCEPTED MANUSCRIPT restitution coefficients, numerical grid resolution, and discretization scheme for convective terms that will always produce a good agreement between numerical and experimental results. Only the 3D simulations allow obtaining physically realistic results that are less sensitive to a choice of numerical parameters.

Keywords: Fluidized beds, validation, 2D, 3D, turbulence, restitution coefficient,

IP

T

specularity coefficient

CR

Nomenclature

US

drag function

CD

d s , d p diameter of solid particle

activation energy

ess

coefficient of restitution

g

gravitational acceleration

g 0,ss

radial distribution function

K

specularity coefficient

kg

turbulence quantities of gas phase

ks

turbulence quantities of solid phase

K pq

interphase momentum exchange coefficient

K gs

interphase momentum exchange coefficient between gas and solid phase

CE

PT

ED

M

AN

EA

nt

AC

k   K-Epsilon turbulence model temperature exponent

p

pressure shared by all phases

ps

solid pressure

R pq

interaction force between phases

Res

Reynolds number for solid particles

RNG renormalization group method

ACCEPTED MANUSCRIPT superficial gas velocity

Ug

volume averaged velocity for gas phase

Vq

volume of phase q

Vsb

volume of the spouted bed reactor

vq

velocity of phase q

vg

velocity of gas phase

vs

velocity of solid phase

vq

particle velocity fluctuation

US

CR

IP

T

U

Greek letters

q

volume fraction of q phase

AN

 s ,max maximum packing limit

volume fraction of solid phase

g

volume fraction of gas phase

q

physical density of phase q

g

physical density of phase g

ˆ q

effective density of phase q

q

q th phase stress-strain tensor

q

bulk viscosity of phase q

s

gs

ED

PT

CE

AC



M

s

collisional energy dissipation kinetic energy transfer

s

granular temperature



fluid-particle interaction coefficient

s

solid shear viscosity

s ,col

solid collisional viscosity

s ,kin solid kinetic viscosity

frictional viscosity

g

viscosity of gas phase

s

volume fraction of solid phase

q

granular temperature

t , g

turbulent viscosity

g

turbulence dissipation rate of gas phase

IP

 s , fr

T

ACCEPTED MANUSCRIPT

k ,g ,  z,g

influence of dispersed phases on continuous phase

standard error of the estimate

US

 est

CR

 k ,  z adjustable constants

p

particle

s

solid phase

M

gas phase

ED

g

Introduction

PT

1

AN

Subscripts

Fluidized beds is a widely used technology in high-temperature conversion of materials including

CE

gasification, waste treatment, and CO2 utilization [1]. Analysis of different types of fluidized beds and their application in industry shows that the next challenging task for engineers is to

AC

increase the efficiency of fluidized beds. This can be done by means of computational fluid dynamics (CFD) modeling of processes occurring inside the fluidized beds. Recently, computer modeling became one of the most crucial elements in the design and optimization of novel technologies in the field of industrial engineering. One of the advantages of computer mode ling is that the behavior and characteristics of the fluidized bed reactor can be investigated without fabrication of a pilot-scale experimental prototype. Therefore, the total costs of product development and design optimization are reduced significantly. However, the numerical simulations of multiphase flows in fluidized beds often involve complex physical and chemical phenomena which have to be modeled using the advanced mathematical models implemented into

ACCEPTED MANUSCRIPT a CFD software. To guarantee adequacy and reliability of any model or software, verification and validation processes are required [2, 3, 4]. The verified and validated numerical models can be used to study the effects of different parameters on the efficiency of the process when experiments are challenging and expensive. With an impressive development in commercial CFD software, e.g., ANSYS FLUENT [5], ANSYS CFX [6], STAR-CCM+ [7], FLOW-3D [8], and in multi-phase flow models, it

T

became possible to simulate the three-dimensional (3D) fluidized beds without oversimplification.

IP

The Eulerian-Eulerian (E-E) and Eulerian- Lagrangian (E-L) are the two main numerical

CR

approaches used in CFD-based fluidized bed models. In the E-E approach, the fluid and solid phases are considered as interpenetrating continuum phases [9]. The conservation equations are

US

required to describe the transportation of mass, momentum, and energy. The E-L approach treats the fluid phase as continuum and the solid as individual dispersed particles. The Newton’s law of

AN

motion and Navier-Stokes equations are used to solve the discrete phase and the fluid phase, respectively [9].

M

Before we proceed with the literature review of the E-E simulations of fluidized beds, we will outline and define the most important parameters (except properties of particles, e.g. density,

ED

size and shape) that determine the hydrodynamics of the gas-solid flow in a fluidized beds, and, therefore, should be accounted for in the numerical model. These parameters are the superficial

PT

gas velocity, the gas-solid drag, the specularity and restitution coefficients, and the frictional viscosity.

CE

The superficial gas velocity represents the inlet velocity of gas through a packed bed in the gasifier. The fluid regime can be changed with the increase of the upward superficial gas velocity

AC

from zero to the point of fluidization. In this case, the upward gas-solid drag force exceeds the downward body force of gravity. The gas-solid drag force (or just the drag) is the resistance force caused by the motion of a particle through a fluid [5]. The drag is the main reason for the bed expansion and intensified mixing between particles and gas. Due to the complex interactions between particle-particle and particle-wall, the description of dynamic behavior of the fluidized bed is complicated. The specularity coefficient is used to characterize particle-wall interaction. The restitution coefficient is used to characterize the particle-particle interaction. Both of them play an important role in the bed behavior. Different hydrodynamic behavior (e.g. flow structure, pressure drop, or particle velocity) can be captured by the model using different values of

ACCEPTED MANUSCRIPT specularity coefficient and restitution coefficient [10, 11]. After collision between the particles, the stress is generated due to the enduring contact and momentum transfer through friction between the particles which is called as frictional viscosity [12]. The CFD-based models to simulate fluidized beds can reveal the influence of each of these parameters and guide the improvement of the process. Even though there are numerous CFD studies exploring the influence of the restitution

IP

T

coefficient ( ess ) and the specularity coefficient ( K ) on the flow behaviour inside the fluidized beds, e.g. see [10, 11], a clear graphical explanation of these coefficients is still needed. For that

CR

reason, we introduce Fig. 1 (a), which shows the definition of the restitution coefficient (particle-particle collision) in the fluidized bed. For perfectly inelastic collision, ess  0 , the

US

particle A and the particle B stick and move together with a final speed in the bed after collision. For perfectly elastic collision, ess  1 , the particle A and the particle B bounce away from each

AN

other and move separately with the same magnitude of final speed after collision. The restitution coefficient ess is defined as follows [13]:

M

v1  v2 u1  u2

(1)

ED

ess =

where v1 and v2 are the post-collision velocities of particle A and particle B, respectively. The

PT

u1 and u2 are the pre-collision velocities of particle A and particle B, respectively. Fig. 1 (b) illustrates the meaning of the specularity coefficient that defines particle-wall

CE

collision in the fluidized bed. For perfectly specular collision, K  0 , the particle collides with the wall and then bounces to the central region of bed with the same incidence angle and the reflection

AC

angle ( 1 =  2 ). For perfectly diffuse collusion, K  1 , the particle has different incidence and reflection angles ( 1   2 ) after collision with the wall. With these definitions in mind, we reviewed published studies on numerical simulations of fluidized beds. Many researchers studied the characteristics of gas-solid fluidized bed reactors by applying the E-E and E- L models in the CFD software, ANSYS FLUENT. Taghipour [14] used the E-E approach to study the effect of different drag sub-models (Syamlal-O’Brien [15], Gidaspow [16], and Wen-Yu [17]) when running the two-dimensional (2D) simulations. The comparison of the numerical and experimental bed pressure drop, bed expansion ratio, and

ACCEPTED MANUSCRIPT gas-flow pattern revealed a good agreement between the numerical and experimental results. Loha et al. [10] carried out an E-E simulation to investigate the effect of specularity coefficient, K , on the hydrodynamic behavior of the bed for 2D simulations. Different values of restitution coefficient ess = (0, 0.01, 0.1, 0.3, 0.6 and 1.0) were tested. It was shown that with the decrease of K from 1.0 to 0.01 , the pressure drop across the bed is constant and remains in

T

the range 1.90 1.92 kPa. However, the bed pressure drop reduces significantly to 1.84 kPa

IP

when K = 0 . The comparison of the numerical and experimental data of time-averaged axial particle velocity and granular temperature showed a good agreement if specularity coefficient is

CR

equal to 0.1 , 0.3 , 0.6 and 1.0 .

Loha [11] used the E-E approach to study the influence of the restitution coefficient, ess ,

US

on the hydrodynamics of the fluidized beds using 2D laminar model. The effect of the restitution coefficient ess in the range (0.85  1.0) on the time-averaged solid volume fraction is not clear.

AN

The assumption of ess  1 results in the flow pattern without any bubbles in the bed which is

M

physically non-realistic. When ess decreases from 1.0 to 0.85 , the quantity and the size of bubbles increases and the averaged solid volume fraction decreases accordingly. The large

ED

fluctuation in the pressure drop within the bed is contributed to the formation of a large number of big bubbles. The numerical results of time-averaged axial particle velocity and granular

PT

temperature were in line with the experimental measurements for ess = 0.95 and 0.99 . Both E-E and E-L models were used by Almohammed [18] to study the effects of the flow

CE

rate of the injected gas on the 3D    dispersed turbulence model. The numerical result showed that the restitution coefficient of 0.6 and specularity coefficient of 0.5 with Syamlal-O’Brien

AC

drag function give the best agreement with the experimental data on the bed height and the equivalent bubble diameter over time. The E-E approach can successfully predict the flow pattern at a mass flow rate of 0.005 kg/s, while a significant deviation is produced as the flow rate increases to 0.006 kg/s. Stroh [9] studied and compared three different numerical methods, Euler-Euler (E-E), Euler-Lagrange-stochastic (E-L-S), and Euler- Lagrange-deterministic (E-L-D) on a 3D laminar model. The comparison of numerical and experimental measured bed expansion and equivalent bubble diameter shows that all three approaches can predict the bubble development for the inlet air mass flow rates of 0.005 and 0.006 kg/s. For the E-E approach, the accurate flow pattern can be

ACCEPTED MANUSCRIPT predicted for the smaller mass flow rate (0.005 kg/s) but not for the higher mass flow rate (0.006 kg/s) between the time interval of 0  300 ms. Cornelissen [19] used the E-E approach to investigate different mesh size, time step, and convergence criteria on a 2D laminar model. The comparison of numerical and experimental data on overall bed voidage shows that the model using 60  500 grid gives a small deviation of 4% in comparison of 45  250 , 50  400 , and 70  700 grids. When time step t = 0.001 s is used, the

IP

T

overall bed voidage values predicted numerically agree well with the experimental values. Jin [20] used the E-E approach to investigate the effect of particle-particle restitution

CR

coefficient, the particle-wall restitution coefficient, and the specularity coefficient on a 3D k   dispersed model. The large particle size ( d p  100 m ) leads to a high sensitivity of

US

particle-particle restitution coefficient on gas and solids velocity and solid volume fraction distribution in comparison to a small particle size ( d p  100 m ). The sensitivity of

AN

particle-particle restitution coefficient on the variable particle density was not explained. The particle-wall restitution coefficient has a negligible effect on the time-averaged solid volume

M

fraction, axial particle and gas velocity. The solids volume fraction decreases from 14% to 4% near

ED

the wall when the specularity coefficient was increased from 0 to 0.01. The comparison between numerical and experimental axial pressure gradient shows that the particle-particle restitution coefficient of 0.99, particle-wall restitution coefficient of 0.99, and specularity coefficient of 0 can

PT

give the best result with an averaged deviation from the experiment of

5.6% .

In series of works Ghoniem and co-workers [21, 22] carried out 2D and 3D numerical

CE

simulations of dense solid–gas cylindrical and rectangular fluidized beds using the E-E based model utilizing the MFIX code. The influence of different wall boundary conditions and drag

AC

models and the specularity coefficient were extensively studied numerically. In particular, it was shown that in the case of thin rectangular fluidized beds the effects of the front and the back walls are significant [22]. The influence of the specularity coefficient, which is responsible for the momentum transfer from the particles to the wall, showed that the appropriate value of the specularity coefficient depends strongly on the superficial gas velocity of the bed. However, it should be noted that no turbulence models were used in that simulation [22]. Li and co-workers [23] investigated numerically the influence of the solid-phase wall boundary conditions on the hydrodynamics of a bubbling fluidized bed. Parametric studies of the

ACCEPTED MANUSCRIPT particle–wall restitution coefficient and specularity coefficient were performed to evaluate their impact on the predicted fluidized bed behavior in terms of bed expansion, local voidage, and solid velocity. Both 2D and 3D simulations (using the MFIX code) were conducted and compared with available experimental data. Good agreement was observed. However, it should be noted that no any turbulence models were used in that simulation. Kong et al. [24] evaluated numerically the effect of wall boundary conditions on the

T

hydrodynamics of circulating fluidized beds. The Eulerian-Eulerian 2D CFD-based model

IP

(ANSYS Fluent 13.0) coupled with k   turbulence model formulated for each phase was used.

CR

The results showed that the specularity coefficient has a strong effect on the solid particle distribution near the wall. At the same time, the particle–wall restitution coefficient has less

US

effect on solid particle distribution near the wall than that of the specularity coefficient. Numerous amount of studies on the Eulerian-Eulerian models applied for the fluidized bed

AN

simulations show that these models are extensively used. The major reason for this is that they are less computationally demanding as compared to other more sophisticated modeling approaches

M

where each particle motion is being resolved. The potential applicability of the E-E approach is high since the unlimited number of particles can be modeled with this approach. The disadva ntage

ED

of the E-E approach is that some discrete character of the solid phase cannot be identified [25]. The additional source terms (closure laws) are required in the E-E approach to define the interactions

PT

between the solid and fluid phases. The advantage of the E-L approach is that the high- level details of the solid phase are available since the individual particles are tracked. On the other hand,

CE

significant computational cost (e.g. CPUs and memory resources) is required for the E-L approach to give a statistically meaningful numerical prediction of the real hydrodynamic behavior [25].

AC

The number of particles that can be modeled using the E-L approach is limited in the range from 106 to 109 [26]. The number of particles in each grid cell is also an important parameter. Due to

the validity of the drag force closure relations used in the E-E approach, the particle size must be smaller than the grid cell size: the edge of the cell in each direction should be at least two particle diameters ( x  2d p , y  2d p , and z  2d p ). Table 1 shows the number of particles in each grid cell for several works. As can be seen, several studies do not follow the requirement o f the particle resolution. These studies are highlighted in the table. The numerical results from these works might not replicate the physical behavior. The Syamlal-O’Brien and Gidaspow drag sub- models cannot be applied in this cases due to less number of particles ( < 10 ) in each cell [26].

ACCEPTED MANUSCRIPT In spite of the numerous works on validation of the CFD-based models (E-E or E-L) applied to simulations of fluidized bed, there is a lack of comparison of the results obtained on 2D and 3D simulations specifically for a thin rectangular bed. Often in the literature, the 2D simulations are used to extrapolate theoretical findings for different types of fluidized beds, thin to wide [27]. However, in the case of a thin rectangular fluidized bed, the wa ll effects in the third dimension might be significant. The purpose of this work is to show the significance of the third

T

dimension on the simulation results. The experimental work of Taghipour et al. [14] is used as a

IP

reference. Our goal is to estimate the magnitude of the difference between 2D and 3D simulations

CR

and compare the numerical and experimental results. A similar effect was investigated by Ghoniem and co-workers [22] but only for the laminar flows. We will consider turbulent models as

US

well. Xie et. al. [28, 29] studied numerically two- and three-dimensional computational modeling of the fluidized beds. A good agreement between 2D and 3D simulations was observed for

AN

bubbling regimes with low gas velocity. For slugging and turbulent regimes, the 2D and 3D simulations produced different results. However, we note that the Authors did not carry out any

M

validation against experimental data.

Motivated by a lack of solid comparison of 2D & 3D simulation results to experimental

ED

data, the main objective of the present study is to validate the 2D & 3D Euler- Euler-based models available in a commercial CFD software, against experimental data [14]. It is also necessary to

PT

analyze the effects of the laminar and turbulent k   flow models on numerical results. Together with the influence of a grid resolution, the sensitivity of key model parameters (the specularity

CE

coefficient, the restitution coefficient, and the drag sub-model) will be evaluated. The rest of the paper is organized as follows. Section 2 presents the description of the

AC

problem. The details of reference experiments are given in this section. The computational model is described in Section 3. Section 4 outlines the results of the present study. Finally, the conclusions are drawn in Section 6.

2

Problem description

To evaluate the capabilities of a commercial CFD software to capture the physics of fluidized beds and study the sensitivity of parameters of a fluidized bed, the experiment of Taghipour [14] is numerically replicated. The principal scheme of the experiment and simulation is shown in Fig. 2.

ACCEPTED MANUSCRIPT According to the experimental setup, the Plexiglas column has a height of 1.5 m, a width of 0.28 m, and a thickness of 0.025 m. Spherical glass particles with a diameter of 250  300  m , density

2500 kg / m3 are fluidized with air at ambient conditions. The static bed height is 0.4 m with a solid volume fraction of 0.6. The maximum packing limit for the solid phase is max  0.63 . In the experiment, the overall pressure drop and bed expansion were monitored at a superficial gas

T

velocity of U = 0.38 m/s and 0.46 m/s. The measurements of time average pressure drop were

IP

recorded at 10 Hz for 20 s intervals once the steady-state conditions were achieved. To compare

CR

the numerical results with the experiment, the voidage measurements were conducted at the height

3

US

0.20 m above the distributor plate across the width of the column.

omputational model

AN

On the basis of the Euler-Euler approach, the unsteady laminar and turbulent Reynolds-averaged Navier-Stokes (RANS) formulations for the mass and momentum are solved for gas and solid

M

phase separately. The gas phase is modelled using the standard k - turbulence model with a standard wall function. The solid phase is modeled with the kinetic theory of granular flows

ED

(KTGF) [5] that describes the fluctuations and collisions between the particles. The momentum equation for the gas and solid phases is given by the modified Navier-Stokes equations that include

PT

the inter-phase momentum transfer terms. In particular, the frictional viscosity, the drag coefficient sub- model, the granular temperature, the bulk viscosities, and the granular temperature

CE

are the additional terms that depict the momentum exchange. Both turbulence interaction and turbulent dispersion between the gas and solid phases are estimated by the Simonin & Viollet

AC

theory [30]. The Tchen theory [31] of dispersion of discrete particles is utilized in the dispersed turbulence model to estimate turbulent quantities (dispersion coefficients, correlation functions, and the turbulent kinetic energy) for the particle phase. The Syamlal-O’Brien [15] and Gidaspow [16] drag sub- models are used. The Syamlal-O’Brien drag sub- model is based on the terminal particle velocities in fluidized beds [32]. The Gidaspow drag model is an integration of the Wen-Yu model and the Ergun model that can be applied for the flow where viscous forces are dominant [33]. As the inlet boundary condition, a superficial gas velocity of 0.38 m/s or 0.46 m/s is set. The outlet boundary condition is a fully-developed flow. The boundary conditions for turbulence

ACCEPTED MANUSCRIPT are specified with uniform value, 5% for turbulence intensity and 0.1 for turbulence viscosity ratio. The boundary conditions in terms of specularity coefficients (0.01, 0.1, 0.5, and 0.9) and the particle-particle restitution coefficients (0.90, 0.95, and 0.99) are applied in the simulation and their effects are studied in detail. The simulations were performed using the CFD commercial software, ANSYS FLUENT 16.2. The computational geometry was discretized with three block-structured grids using a finite

T

volume method with 56  300 , 112  600 , and 56  300  5 cells (see Table 2). The particle

IP

diameter considered in simulations was d p = 275  m .

CR

The pressured-based coupled scheme [5] was used to solve mass and momentum conservation equations corresponding to the behavior of the fluidized bed. To fulfill the coupling

US

algorithms, the implicit discretization of both pressure gradient terms and face mass flux were included in the equation. Quadratic upwind interpolation for convection kinematics (QUICK)

AN

scheme [34] was activated to discretize the convection terms in the momentum equations. The modified type of the high-resolution interface capturing (HRIC) scheme [35] was used to estimate

M

the volume fraction of gas or solid. In comparison between the QUICK and the first-order upwind schemes, the accuracy of numerical solution can be improved significantly by the modified HRIC

ED

scheme since the non- linear blend of upwind and downwind differencing are included in this scheme [35]. The second-order upwind scheme for the discretization of the convective terms was

PT

used in the k - model. The details of the models and schemes used in the numerical simulation of fluidized beds are given in Table 3.

CE

The total physical simulation time of 20 seconds was completed for each simulation. The time-averaged values were obtained from the last 18 seconds after a dynamic steady state

AC

converged solution was reached within two seconds after the process initiation. To avoid instability, a small constant time step of 104 s was used for each case, and the number of inner iterations on each time step was set to 50. The total of 31 simulations (eight 3D cases and twenty three 2D cases) were carried out in this study. The cases are described in Table 4. Finally, it should be noted that in all 2D and 3D simulations we used the full bed according to its sizes in 2D and 3D cases.

4

Results

ACCEPTED MANUSCRIPT To gain a general understanding of the accuracy of the numerical predictions of voidage fraction, the standard deviation profiles of the 2D and 3D simulations were estimated and are shown in Fig. 3. The formula of standard deviation [36] is shown below:

 est % =

(Y

simulation

 Yexperiment )2

N

100

(2)

T

where  est is the standard error of the estimate, Yexperiment and Ysimulation are the experimental and

IP

predicted points of the voidage fraction along x  coordinate at z = 0.2 m, respectively. The standard deviation provides reliable measure of error as compared to other functions (e.g. weight

CR

standard error) since this method comes from Linear Least Squares Regression used by curve fitting. The standard deviation of 3D and 2D models is in the range [3.85  8.54] % and

US

[3.97  13.08] %, respectively. The difference between the maximum and minimum deviations for

3D and 2D simulations is 4.69% and 9.10%. This reveals that the results of the 3D simulations of

AN

fluidized bed agree better with the experimental results, and therefore, should be used. However, a proper choice of parameters in 2D simulations might also lead to results that match experiments.

M

The detailed analyses of the results is presented in the following order. First, we assess the

ED

difference between two- and three-dimensional simulations. Second, the results for the laminar and turbulent models are discussed. Effects of the coefficient of restitution and specularity

4.1

CE

viscosity sub-models.

PT

coefficients are considered next. Finally, we examine the effects of the drag force and the frictional

Three-dimensional versus two-dimensional approaches

AC

This section compares the numerical predictions of two different modelling approaches, 3D and 2D. For a 3D model, the total simulation time of 2 weeks is required to simulate the physical process time of 20 seconds. However, a 2D model needs only 2

3 days to simulate the same

physical time. To illustrate the fluidized bed behavior, Fig. 4 shows the snapshots at different time instances of the slice of a 3D contour plot of the solid volume fraction and the turbulent viscosity ratio. The red color represents the maximum solid volume fraction and turbulent viscosity ratio, whereas the blue color denotes the minimum value of zero. As can be seen from Fig. 4 (a), bubbles formed at the bottom of the bed move up to the bed surface with the increasing bubble size. The

ACCEPTED MANUSCRIPT shape of the bubbles ceases to remain spherical because of the intensive breakage and coalescence. Upon reaching the bed surface, the solid particles are forced towards the wall with the burst of big bubbles. Then the particles fall down along the walls due to the gravitational force. This behaviour of particles can be seen in Fig. 5 that presents a time-averaged solid volume fraction in the middle slice of the domain. A lower solid volume fraction (30%) is observed in the central area of the

T

bed while the higher solid volume fraction of 60% is achieved in the regions close to the wall. The time snapshots of the turbulent viscosity ratio for the same 3D simulation case are

IP

depicted in Fig. 4 (b). The turbulent viscosity ratio – defined as a ratio between turbulent t and

CR

molecular o dynamic viscosities – is necessary to estimate the level of turbulence within the simulation domain. As can be seen from the Fig. 4, the distribution of turbulence is consistent with

US

the gas flow patterns. High turbulent viscosity ratio (t / o > 10) exists inside the big bubble space at the startup phase ( t < 2.0 s). After reaching the dynamic steady state ( t > 2.0 s), weak

AN

turbulent viscosity ratio in the range (0  10) is observed inside the bed, while strong turbulence

M

in the range (15  110) is distributed above the bed surface. This indicates that only a small amount of gas (bubbles) is trapped inside the bed. The gas primarily fills the space above the bed

ED

surface by releasing the trapped gas from the bubbles. To study the transient dynamics of the gas and solid flows in a fluidized bed, we use the

flow velocity U s [37]:

PT

volume-averaged or global velocity of the solid phase in the entire bed. The volume-averaged solid

1 V fb

CE

Us =

x

y

 0

0

H

0





us2, x  us2, y  us2,z d x d y d z

(3)

AC

where V fb is the volume of the fluidized bed reactor and us is the velocity of the solid phase. This global velocity is used to calculate a spinup or startup time for the volume- force driven flows [37]. Time histories of U s calculated using different values of specularity coefficient for 2D and 3D RANS approaches are depicted in Fig. 6. A noticeable narrow peak in the volume-averaged velocity profile within the first 2 seconds is captured in all cases. This peak corresponds to the startup time. Therefore, the fluidized bed reaches a dynamic steady state within two seconds of process initiation. After the startup period, steady fluctuations of solid volume fraction are observed in the bed surface that might be explained by continuous formation, rise, coalescence, and the breakup of bubbles.

ACCEPTED MANUSCRIPT Analysis of U s time history shown in Fig. 6 reveals that the solid velocities predicted using the 2D model are higher than velocity obtained from the 3D model. Dominant velocity in the ranges 0.2  0.3 m/s and 0.1  0.2 m/s are found for 2D and 3D cases, respectively. The lower solid velocity from the 3D models is explained by the loss of extra friction due to the increased contact area between the fluid and the walls. To capture this phenomenon, the simulation of the fluidized bed has to use the 3D approach rather than the 2D approach to obtain physically realistic

IP

T

results. Even if the 2D grid is refined from 56  300 to 112  600 , the numerical results do not approach the experimental measurements. However, significant computational demands

CR

associated with the 3D simulations make 2D simulations more appealing for the sensitivity analyses of the parameters. For that reason, it is necessary to assess to what extent the 2D results

US

can be used for these purposes.

Since the main goal of this study is to validate the numerical approach, we will start

AN

analyses of the 3D numerical results that replicate the real phenomena better than 2D simulations. Fig. 7(a) compares the predicted void fraction values against the experiment data along x 

M

coordinate at z = 0.2 m with the inlet gas velocity of 0.38 m/s. The experimental data is adapted from the experiment carried out by Taghipour [14]. The averaged deviation from the experimental

ED

values for all three numerical calculations in Fig. 7 is 5.36% suggesting that the 3D approach produces a good agreement with the experimentally measured void fraction profiles. It can be seen

PT

that there is a flat distribution with high gas volume fraction in the central region and sharp distribution with low gas volume fraction close to the wall. The axial variations in gas volume

CE

fraction are similar for the coefficient of restitution 0.9 and 0.95 without activation of the fictional viscosity sub- model. The deviations from experimental data are 4.86% and 4.71% for these cases,

AC

respectively. It confirms that there is no qualitative change in the nature of variation of the results by varying the coefficient of restitution within a range of [0.9  0.95] . However, a relatively large deviation of 6.50% against experimental data is presented by the increase of specularity coefficient, K , from 0.1 to 1. This overestimated prediction is caused by a perfectly diffuse collision between the particle and the walls when K = 1 . Fig. 7 (b) depicts the 3D simulated time-averaged void fraction profile with the inlet gas velocity of 0.46 m/s. By comparing the numerical results with the experiment points, all four cases give small deviation in the range between 3.85% to 6.37%. With the increase of the specularity coefficient from 0.1 to 1, the minimum value of void fraction (41.5%) can be observed on the wall, the main possible reason is

ACCEPTED MANUSCRIPT explained in the previous discussion. Finally, Fig. 7 (c) presents time-averaged void fraction profiles predicted using 3D laminar model and 3D k - disperse turbulence model. It can be seen that the turbulence model produced better agreement with experiments.

4.2

Laminar model versus RANS model

T

The three-dimensional simulations of fluidized bed using an unsteady RANS k   model show

IP

that the flow within the bed is laminar which is indicated by a low turbulent viscosity ratio (< 10)

CR

. At the same time, stronger turbulence develops above the bed surface. Therefore, it is necessary to assess what flow model (laminar or turbulent) should be used. The 2D fluidized bed configuration was simulated on grids Grid1 56  300 and Grid2 112  600 . The particle

US

concentration distributions inside the fluidized bed for the laminar and RANS models are shown in Fig. 8 for Grid1 and in Fig. 9 for Grid2, respectively. Similar to the results of the 3D simulation, in

AN

2D simulations, the bubbles undergo the formation, rise, coalescence, and burst. By comparing the flow pattern in both laminar (see Fig. 8 (a) and 9 (a) and RANS (see Fig. 8 (b) and 9 (b)) models, it

M

is noticed that the solid flow patterns inside the bed are relatively symmetrical for both models at

ED

the startup state (t < 2 s) . After reaching the dynamic steady state, similar bubble size and expanded bed height can be observed for both laminar and RANS models. In order to study the overall behavior of the fluidized bed, the 2D time-averaged solid

PT

fraction contour plots and gas velocity vector plots for laminar and RANS models on different grids are depicted in Fig 10 and 11. The symmetry in the solid pattern and gas velocity is achieved

CE

in all cases. The solid-phase distribution is similar for both laminar and RANS models in both girds suggesting that fluidization hydrodynamics is not sensitive to the choice of the flow model. A

AC

higher solid volume fraction of 0.6 is achieved close to the wall and at the bottom of the bed (0.45). A lower solid volume fraction of 0.3 is observed in the center of the bed. This distribution of solid volume fraction reveals that more bubbles are formed at the center region rather than at the bottom. The 15% solid volume fraction is observed at the bed surface demonstrating the burst of bubbles when they approach the surface. The analyses o f the distribution of fluidizing gas velocity (Fig. 10 (b) and Fig. 10 (b)) reveal that the fluidizing gas velocity in the center of the bed is high and gradually decreases above the bed surface. When the fluidizing gas flows down to the bottom, it changes direction after encountering the inlet gas and the magnitude of the velocity vector

ACCEPTED MANUSCRIPT decreases. Consequently, the mixture of fluidizing gas and inlet gas forms a vortex flow of gas in both sides between the center and the wall. This vortex is an important feature of the fluidized beds due to its governing role in heat and mass transfer applications. In particular, the residence time of solid particles in this zone is significantly larger than in the bulk zone. Fig. 12 displays the contour plots of the turbulent viscosity ratio predicted at different time instances using the 2D geometry with 56  300 and 112  600 mesh resolution. The observations

T

are similar to the 3D results: weak turbulence distributed in the bed, and high turbulence exists

IP

above the bed surface. In the startup state, we barely can see turbulence in the entire domain. After

CR

reaching the dynamic steady state, the turbulent viscosity ratio in the packed bed is in the region between 0 and 52, while in the gas phase is in between 103 and 310 for the model. As in the 3D

US

cases, the possible explanation is that less gas is being present in the bed. To illustrate qualitative differences between the laminar model and RANS model, Fig. 13

AN

depicts the time-averaged velocity of the gas phase along the center line corresponding to

x = 0.15 m. Higher gas velocity ( 0.25  1.7 m/s) exits in the bed while lower gas velocity (

M

0  0.5 m/s) distributed in the gas phase for both laminar model and RANS model. This can be explained by the presence of particles. In the gas phase, we obtain lower gas velocity between 0

ED

m/s to 0.38 m/s in the 2D RANS model and relatively higher velocity between 0.15 m/s to 0.44 m/s in the 3D RANS model. The reason for this phenomenon is that higher turbulent viscosity ratio

PT

(112  450) exists in the 2D model and the lower turbulent viscosity ratio (37-110) is distributed

in the 3D model.

CE

Referring to the comparison between solid phase velocities calculated using the laminar and turbulent models, Fig. 14 (a), (b) shows the time histories of the volume-averaged velocity of

AC

the entire solid phase predicted with laminar and RANS models, respectively. It can be seen that both laminar and RANS models show a similar trend of the solid velocity in the 2D model (in the range 0.2  0.3 m/s) and 3D model (in the range 0.15  0.25 m/s), which reveals that the turbulence does not play a significant role on solid velocity. Based on our simulation results that compare the laminar and turbulent models applied for the fluidized bed simulations, we conclude that the laminar model can be used when the isothermal fluidized bed is considered. For gasification and combustion applications where chemical reactions take place, the RANS model should be used. In these cases, species transport is affected by intensive turbulence above the bed surface since diffusion is different as compared to the

ACCEPTED MANUSCRIPT laminar flow.

4.3

Effect of coefficient of restitution ess

Fig. 15 (a) presents the axial distributions of the time-averaged void fraction for different values of coefficient of restitution ( ess = 0.90 , 0.95 , and 0.99 ) obtained as result of the 2D simulations.

T

The core-annular structure of the flow is predicted when collisions are inelastic. The void fraction

IP

is larger in the central region and less in the near wall region. A higher void fraction with a value of

CR

46% is observed on the wall for ess = 0.90 . This can be explained by the formation of many larger bubbles at a lower value of coefficient of restitution. The deviation of numerical results against

US

experiment data are 6.17%, 7.60%, and 7.72% for ess = 0.90,0.95, and 0.99 , respectively. This

4.4

AN

suggests that the void fraction is not sensitive to the choice of ess in the range [0.90  0.99] .

Effect of specularity coefficient K

M

Four different values of the specularity coefficients ( K ) in the range [0.01  0.9] are considered

ED

in the 2D numerical models. The comparison of the experimental and numerical data is presented in Fig. 15 (b). The differences in the numerical results using various specularity coefficients are more distinct while studying the time-averaged void fraction distribution. With the decrease of the

PT

specularity coefficient from 0.9 to 0.1, the trend of the void fraction with the deviation in the range

CE

between 7.57% and 7.99% can be observed. This suggests that the influence of K in the range between 0.1 to 0.9 on the numerical prediction of void fraction pattern is negligible. The existence of a core-annular structure of the flow was predicted for K between 0.1 and 0.9. When the

AC

specularity coefficient decreases to 0.01, a different gas volume fraction pattern is observed with the downflow of gas near the wall and at the center with the deviation of 9.68%. Two maximum gas volume fraction values, around 0.64, can be found at x = 0.05 and 0.225 m. The lower value of the specularity coefficient means the resistance between the particle and the wall is small. Consequently, more concentrated particles are distributed close to the wall, resulting in two peaks of the void fraction located each side between the center and the wall. Fig. 16 reveals that the solid velocity is significantly affected by the variation of specularity coefficient in the range between 0.1 and 1.0. For the superficial gas velocity of 0.38 m/s, the solid

ACCEPTED MANUSCRIPT velocity pattern of K = 1 model is lower than that of K = 0.1 model. The dominant solid velocity is around 0.1 m/s and 0.15 m/s for each model. The unity value of specularity coefficient leads to zero fluid velocity on the wall resulting in such a lower value of the time-averaged solid velocity. The effect of superficial gas velocity is more pronounced for the solid velocity profile. With the increase of the superficial gas velocity from 0.38 m/s to 0.46 m/s, the higher solid velocity pattern with an averaged value of 0.2 m/s is observed. Particles are carried up in the

T

central region and fall down in the near wall region by the inlet gas. Thus, a positive relation

4.5

CR

IP

between the superficial gas velocity and the particle velocity results in this phenomenon.

Effect of drag force and frictional viscosity

US

Fig. 17 (a) presents a comparison of the predicted time-averaged gas volume fraction with the experimental results by different parameters using Grid1. A better numerical result with a

AN

deviation of 4.50% is contributed from the RANS model by taking into account the Symal-Obrien drag sub- model and Johnson-et-al frictional viscosity. Whereas the laminar model gives the largest

M

discrepancy, 9.30%, with the experimental measurements. The simulation using Syamlal-Obrien drag submodel gives 6.17% deviation, while that using Gidaspow sub- model shows a deviation of

ED

8.78%. This reveals that use of the Syamlal-Obrien drag sub-model to predict the gas flow pattern in this fluidized bed gives more accurate results. The RANS model with a Syamlal-Obrien drag

PT

sub- model presents a relatively fully-developed bubbly flow pattern while the laminar model shows an annular flow pattern from the profile. This can be explained by irregular fluid flow from

CE

the turbulent regime and uncrossed orderly fluid flow produced by the laminar regime. The velocity of the fluid is not constant at any point of fluid for the turbulent regime, but constant

AC

velocity can be represented by the laminar regime. Consequently, the distribution of the void fraction pattern is consistent with the fluid velocity profile. Fig. 17 (b) shows the numerical results based on Grid2. The profile depicts that the simulation with the frictional viscosity can provide better agreement with the experimental data. The deviation for the frictional viscosity sub- model is 6.19%, while the deviation for the RANS model without applying frictional viscosity is 10.15%. A similar trend of the time-averaged void fraction profile can be predicted by different models without taking into account the frictional viscosity. The RANS model with Syamlal-Obrien drag submodel is 1.64% accurate than in comparison to both the Laminar model and RANS model with Gidaspow drag sub-model. The

ACCEPTED MANUSCRIPT largest difference between predicted gas volume fractions and experimental measurement appeared in the distance interval between 0.15 and 0.28 meters. This can be attributed to the fact that the choice of appropriate multiphase model or the drag sub-model along might probably not improve the numerical results. Some other parameters such as frictional viscosity sub- model, specularity coefficient, and wall roughness are not properly modeled in the bed. The influence of frictional viscosity is depicted using both grids. The RANS model without

T

frictional viscosity demonstrates an averaged excessive prediction, 8.16% across the full

IP

cross-section, whereas the model with frictional viscosity shows an averaged deviation of 5.35%.

CR

The stress that is generated by friction between the particles can be contributed to the solids shear viscosity. Consequently, higher energy loss and lower gas volume fraction appeared due to the

Discussion

AN

5

US

high viscous solid particles.

Before moving to the conclusion of the present study, it is important to discuss the appropriateness

M

of 2D simulations to validate different numerical models and sub- models. The number of

ED

publications reporting the results of the 2D modeling of fluidized beds, e.g. see works [27, 19, 38, 39] and a list of references within these publications, tremendously overweights a limited number of the 3D studies. The reason is well understood: the 3D simulations are associated with significant

PT

computational demands and time. Xie et. al. [28, 29] showed that in some cases the similarities between 2D and 3D simulations can be observed for bubbling regimes with low gas velocity. Hu

CE

et. al. [39] used the 2D simulations to validate their EMMS drag model. Excellent agreement with the experimental data was demonstrated. However, analysis of the results of our present work

AC

shows that the 2D model does not produce physically meaningful results in terms of void fraction distribution. This discrepancy is attributed to neglecting the third dimension in 2D simulations. The results of the 3D simulations clearly show that wall friction in the third direction plays a significant role in the size of bubbles and in the spatial distribution of the volume fraction of particles inside the fluidization zone. We show that a reasonable set up of model parameters (without any sophisticated new drag models) produces very good agreement with the experimental data when 3D simulations are used. We emphasize that a validation of new sub- models against experimental data using 2D simulations does not guarantee that such sub- models work well in 3D

ACCEPTED MANUSCRIPT simulations giving physically realistic results. Moreover, a certain combination of numerical parameters for the 2D simulations can result in acceptable agreement with the experimental data which could be a misleading finding. We suggest that 3D models are used to validate the newly-developed sub- models and to select the model parameters. Only then the validated sub-models and an outlined set of parameters can be used for 2D sensitivity analyses. Finally, we want to discuss the role of turbulence and its modeling on the results presented

T

in this work. We used a well accepted and well used k - dispersed turbulence model, which

IP

showed a reasonable agreement with the experimental studies [18]. The analysis of literature

CR

carried out in section Introduction showed that many 2D and 3D CFD simulations of fluidized beds using laminar flow model, e.g., see [23, 21, 22], showed reasonable agreement with the

US

experimental data published in the literature. Our simulations revealed that the turbulence inside the bed is significantly weaker than in the gas phase outside of the bed. However, there are

Conclusions

M

6

AN

numerous works devoted to fundamental studies of the turbulence in fluidized beds [40, 41].

ED

In the present study, the numerical simulations of a fluidized bed filled with the Geldart Group B particles were performed using the CFD commercial software, FLUENT 16.2. The two- and three-dimensional approaches were considered. A gas-solid phase Eulerian-Eulerian based model

PT

was used to replicate the experimental of Taghipour[14]. The influence of the model dimension (2D & 3D), the flow regimes model (laminar & turbulent), the model parameters (the restitution

CE

coefficient responsible for particle-particle interaction and the specularity coefficient responsible for the wall-particle interaction) and the computational girds were studied. Following are the main 

AC

conclusions that can be drawn from this study: The 3D and 2D approaches showed a deviation of a time-averaged void fraction within [3.85-8.54]% and [3.97-13.08]%, respectively. 

The results of the 3D simulations clearly illustrate that the finite size of the domain where fluidization occurs plays a significant role in the size of bubbles and in spatial distribution of the volume fraction of particles inside the fluidization zone.



We showed that the 2D approach does not give a good agreement with the experimental data. However, some combination of the grid resolution, the specularity and restitution

ACCEPTED MANUSCRIPT coefficients can give results close to the experimental points. 

It was shown that the turbulence does not play a significant role inside the fluidized bed. However, the use of    model provides slightly better agreement with experimental data in comparison with results given by a laminar model.



Taking into account the frictional viscosity provides better agreement with the

The specularity coefficient has a significant influence on the voidage distribution inside the

IP



T

experimental data.

bed in comparison to the restitution coefficient under conditions used in the experiment by

CR

Taghipour[14].

Christopher Higman and Maarten Van der Burgt. Gasification. Gulf professional

AN

[1]

US

References publishing, 2011. [2]

P. Roache, K. Ghia, and F. White. Editorial policy statement on the control of numerical

P. Roache. Verification and Validation in Computational Science and Engineering.

ED

[3]

M

accuracy. ASME J. Fluid Eng, 108:1, 1986.

Hermosa Publishers, Albuquerque, NM, 1998. P. Roache. Perspective: Validation{what does it mean? ASME J. Fluid Eng, 131:034503,

PT

[4]

2008.

Ansys Fluent. 12.0 theory guide. Ansys Inc, 5, 2009.

[6]

CFX Ansys. Solver theory guide. Ansys CFX Release, 11:1996{2006, 2006.

[7]

USER CD-adapco. Guide: Star-ccm+, 2009.

[8]

CW Hirt and B Nichols. Flow-3D user's manual. Flow Science Inc, 1988.

[9]

Alexander Stroh, Falah Alobaid, Max Thomas Hasenzahl, Jochen Hilz, Jochen Ströhle,

AC

CE

[5]

and Bernd Epple. Comparison of three different CFD methods for dense fluidized beds and validation by a cold flow experiment. Particuology, 2016. [10]

Chanchal Loha, Himadri Chattopadhyay, and Pradip K Chatterjee. Euler-Euler CFD modeling of fluidized bed: Influence of specularity coefficient on hydrodynamic behavior. Particuology, 11(6):673–680, 2013.

ACCEPTED MANUSCRIPT [11]

Chanchal Loha, Himadri Chattopadhyay, and Pradip K Chatterjee. Effect of coefficient of restitution in euler–euler CFD simulation of fluidized-bed hydrodynamics. Particuology, 15: 170–177, 2014.

[12]

David G Schaeffer. Instability in the evolution equations describing incompressible granular flow. Journal of differential equations, 66(1):19–50, 1987. Peter McGinnis. Biomechanics of sport and exercise. Human Kinetics, 2013.

[14]

Fariborz Taghipour, Naoko Ellis, and Clayton Wong. Experimental and computational

T

[13]

IP

study of gas–solid fluidized bed hydrodynamics. Chemical Engineering Science,

CR

60(24):6857–6867, 2005.

Madhava Syamlal and Thomas J O'Brien. Computer simulation of bubbles in a fluidized

US

[15]

bed. In AIChE Symp. Ser, volume 85, pages 22–31, 1989. Dimitri Gidaspow. Multiphase flow and fluidization: continuum and kinetic theory descriptions. Academic press, 1994.

CY Wen and Y Yu. Mechanics of fluidization. In Chem. Eng. Prog. Symp. Ser., volume 62,

M

[17]

page 100, 2013.

Naser Almohammed, Falah Alobaid, Michael Breuer, and Bernd Epple. A comparative

ED

[18]

AN

[16]

study on the influence of the gas flow rate on the hydrodynamics of a gas–solid spouted

PT

fluidized bed using Euler–Euler and Euler–Lagrange/DEM models. Powder Technology, 264:343–364, 2014.

J.T. Cornelissen, F. Taghipour, R. Escudié, N. Ellis, and J.R. Grace. CFD modelling of a

CE

[19]

liquid–solid fluidized bed. Chemical Engineering Science, 62(22):6334–6348, 2007. Baosheng Jin, Xiaofang Wang, Wenqi Zhong, He Tao, Bing Ren, and Rui Xiao. Modeling

AC

[20]

on high-flux circulating fluidized bed with Geldart group B particles by kinetic theory of granular flow. Energy & fuels, 24:3159–3172, 2010. [21]

A Bakshi, C Altantzis, RB Bates, and AF Ghoniem. Eulerian–Eulerian simulation of dense solid–gas cylindrical fluidized beds: Impact of wall boundary condition and drag model on fluidization. Powder Technology, 277:47–62, 2015.

[22]

C. Altantzis, R.B. Bates, and A.F. Ghoniem. 3D Eulerian modeling of thin rectangular gas-solid fluidized beds: Estimation of the specularity coefficient and its effects on bubbling dynamics and circulation times. Powder Technology, 270:256–270, 2015.

ACCEPTED MANUSCRIPT [23]

Tingwen Li, John Grace, and Xiaotao Bi. Study of wall boundary condition in numerical simulations of bubbling fluidized beds. Powder Technology, 203:447–457, 2010.

[24]

Lei Kong, Chao Zhang, and Jesse Zhu. Evaluation of the effect of wall boundary conditions on numericalsimulations of circulating fluidized beds. Particuology, 13:114–123, 2014.

[25]

Wenqi Zhong, Aibing Yu, Guanwen Zhou, Jun Xie, and Hao Zhang. CFD simulation of

T

dense particulate reaction system: Approaches, recent advances and applications.

Martin Anton van der Hoef, M van Sint Annaland, NG Deen, and JAM Kuipers. Numerical

CR

[26]

IP

Chemical Engineering Science, 140:16–43, 2016.

simulation of dense gas-solid fluidized beds: a multiscale modeling strategy. Annu. Rev.

[27]

US

Fluid Mech., 40:47–70, 2008.

N. Yang, W. Wang, W. Ge, and J. Li. CFD simulation of concurrent-up gas-solid flow in

AN

circulating fluidized beds with structure-dependent drag coefficient. Chemical Engineering Journal, 96:71–80, 2003.

N. Xie, F. Battaglia, and S. Pannala. Effects of using two- versus three-dimensional

M

[28]

computational modeling of fluidized beds part i, hydrodynamics. Powder Technology,

ED

182:1–13, 2008.

N. Xie, F. Battaglia, and S. Pannala. Effects of using two- versus three-dimensional

PT

[29]

computational modeling of fluidized beds: Part ii, budget analysis. Powder Technology,

[30]

CE

182:14–24, 2008.

C Simonin and PL Viollet. Predictions of an oxygen droplet pulverization in a

AC

compressible subsonic coflowing hydrogen flow. Numerical Methods for Multiphase Flows, FED91, pages 65–82, 1990. [31]

JO Hinze. Turbulence, mcgraw-hill publishing co. New York, NY, 1975.

[32]

Sebastian Zimmermann and Fariborz Taghipour. CFD modeling of the hydrodynamics and reaction kinetics of fcc fluidized-bed reactors. Industrial & engineering chemistry research, 44 (26):9818–9827, 2005.

[33]

Ling Zhou, Lingjie Zhang, Ling Bai, Weidong Shi, Wei Li, Chuan Wang, and Ramesh Agarwal. Experimental study and transient CFD/DEM simulation in a fluidized bed based on different drag models. RSC Advances, 7(21):12764–12774, 2017.

ACCEPTED MANUSCRIPT [34]

Brian P Leonard. A stable and accurate convective modelling procedure based on quadratic upstream interpolation. Computer methods in applied mechanics and engineering, 19(1):59–98, 1979.

[35]

S Muzaferija, M Peric, P Sames, and T Schellin. A two-fluid navier-stokes solver to simulate water entry. In Proceedings of the 22nd symposium on naval hydrodynamics, Washington, DC, pages 277{289, 1998. Henry Scheffe. The relation of control charts to analysis of variance and chi-square tests.

T

[36]

P. Nikrityuk, M. Ungarish, K. Eckert, and R. Grundmann. Spin-up of a liquid metal flow

CR

[37]

IP

Journal of the American Statistical Association, 42(239):425–431, 1947.

driven by a rotating magnetic field in a finite cylinder: A numerical and an analytical study.

[38]

US

Phys. Fluids, 17:067101, 2005.

Mehran Ehsani, Salman Movahedirad, and Shahrokh Shahhosseini. The effect of particle

AN

properties on the heat transfer characteristics of a liquid–solid fluidized bed heat exchanger. International Journal of Thermal Sciences, 102:111–121, 2016. S. Hu, X. Liu, N. Zhang, J. Li, W. Ge, and W. Wang. Quantifying cluster dynamics to

M

[39]

improve EMMS drag law and radial heterogeneity description in coupling with gas-solid

[40]

ED

two-fluid method. Chemical Engineering Journal, 307:326–338, 2017. Jesse Capecelatro, Olivier Desjardins, and Rodney O. Fox. Strongly coupled fluid-particle

PT

flows in vertical channels. i. reynoldsaveraged two-phase turbulence statistics. Physics of Fluids, 28: 033306, 2016.

S. Schneiderbauer. Validation study on spatially averaged two-fluid model for gas-solid

CE

[41]

flows: Ii. application to risers and fluidized beds. AIChE Journal, 64:1606–1617, 2018. Adnan Almuttahar and Fariborz Taghipour. Computational fluid dynamics of a circulating

AC

[42]

fluidized bed under various fluidization conditions. Chemical Engineering Science, 63(6):1696–1709, 2008. [43]

Leina Hua, Junwu Wang, and Jinghai Li. CFD simulation of solids residence time distribution in a cfb riser. Chemical Engineering Science, 117:264–282, 2014.

[44]

Tingwen Li, Aytekin Gel, Sreekanth Pannala, Mehrdad Shahnam, and Madhava Syamlal. CFD simulations of circulating fluidized bed risers, part i: Grid study. Powder Technology, 265:2–12, 2014.

ACCEPTED MANUSCRIPT [45]

Thanh DB Nguyen, Myung Won Seo, Young-Il Lim, Byung-Ho Song, and Sang-Done Kim. CFD simulation with experiments in a dual circulating fluidized bed gasifier. Computers & Chemical Engineering, 36:48–56, 2012.

[46]

Regis Andreux, Geoffrey Petit, Mehrdji Hemati, and Olivier Simonin. Hydrodynamic and solid residence time distribution in a circulating fluidized bed: experimental and 3D computational study. Chemical Engineering and Processing: Process Intensification,

CKK Lun, S Br Savage, DJ Jeffrey, and N Chepurniy. Kinetic theories for granular flow:

IP

[47]

T

47(3):463–473, 2008.

Journal of fluid mechanics, 140:223–256, 1984.

Paul C Johnson and Roy Jackson. Frictional–collisional constitutive relations for granular

US

[48]

CR

inelastic particles in couette flow and slightly inelastic particles in a general flowfield.

materials, with application to plane shearing. Journal of fluid Mechanics, 176:67–93,

[49]

AN

1987.

Jianmin Ding and Dimitri Gidaspow. A bubbling fluidization model using kinetic theory of

M

granular flow. AIChE journal, 36(4):523–538, 1990.

ED

Table 1: Summary of published studies. The ratios x / d p , y / d p , and z / d p represent the number of particles along each grid edge; E-E = Eulerian-Eulerian; E-L = Eulerian-Lagrangian; E-L-S = Eulerian-Lagrangian-stochastic; E-L-D = Eulerian-Lagrangian-deterministic x / d p , y / d p , z / d p

E-E, E-L-S, E-L-D

0.75, 0.36, 0.93

Gray Almohammed [18]

E-E, E-L

0.75, 0.36, 1.85

Gray Ehsani [38]

E-E

0.33, 0.33

Cornelissen [19]

E-E

1.87, 1.77

Almuttahar [42]

E-E

14.5, 283

Loha[11]

E-E

3.66, 24.3

Bakshi [21]

E-E

25.1, 17.3

AC

CE

Gray Stroh [9]

PT

Model

Authors

8.30, 5.0 20.6, 17.3 Hua [43]

E-E

65.5, 65.5, 453

Li [44]

E-E

45.7, 94.1

ACCEPTED MANUSCRIPT 11.4, 23.5 22.8, 47.1 20, 20 14.1, 14.1 10, 9.98 E-E

20, 20

Andreux [46]

E-E

157.1, 412.1, 157.1

IP

T

Nguyen [45]

CR

Table 2: List of different grid resolutions used in simulations of fluidized beds. Number of cells

x / d p , y / d p , z / d p

D-Grid1

56  300 = 16800

18.18, 18.18

D-Grid2

112  600 = 67200

D

56  300  5 = 84000 18.18, 18.18, 18.18

US

Grid

AN

9.09, 9.09

M

Table 3: List of different models and schemes used in simulations of fluidized beds. Name

Model/ Scheme Name Eulerian-Eulerian [5]

ED

Multiphase Flow

Implicit Scheme [5]

Viscous Model

Laminar & RANS- k - model (dispersed) [5]

Drag Model

PT

Volume Fraction Parameters

Granular Temperature

Syamlal-O’Brien & Gidaspow [15], [16], [5] Algebraic [5] Syamlal-O’Brien [15], [5]

Granular Bulk Viscosity

Lun et al. [47], [5]

Frictional Viscosity

Johnson-et-al [48], [5]

Frictional Pressure (pascal)

Based-ktgf [49], [5]

Solid Pressure

Syamlal-O’Brien [15], [5]

Radial Distribution

Syamlal-O’Brien [15], [5]

Pressure-Velocity Coupling

Coupled Scheme [5]

Spatial Discretization-Gradient

Least Squares Cell Based [5]

Spatial Discretization-Momentum

QUICK [34]

Spatial Discretization-Volume

Modified HRIC [35]

AC

CE

Granular Viscosity

Fraction

ACCEPTED MANUSCRIPT k - 2nd order Upwind [5]

Transient Formulation

First Order Implicit

Flow

velocity U,

model

Drag force

m/s 3D-Grid

0.38

Laminar

1 0.38

3D-Grid

0.38

1 4

3D-Grid

0.38

1 5

3D-Grid

0.46

1 3D-Grid

0.46

0.1

Syamlal-O’

k -

Brien

RANS

Syamlal-O’

k -

Brien

RANS

Syamlal-O’

k -

Brien

RANS

Syamlal-O’

k -

Brien

0.9

0.95

Syamlal-O’

0.9

0.9

Johnson-et-a

0.1

None

0.1

None

1

Johnson-et-a

l

l 0.1

Johnson-et-a l

0.9

0.5

None

0.9

0.1

None

0.9

1

Johnson-et-a

Brien

RANS

Syamlal-O’

k -

Brien

2D-Grid 0.38

RANS

Gidaspow

0.9

0.1

None

1

k -

0.38

Laminar(H

Syamlal-O’

0.9

0.1

None

RS)

Brien

Laminar(H

Syamlal-O’

0.9

0.1

Johnson-et-a

RS)

Brien

Laminar(U

Syamlal-O’

DS)

Brien

RANS

Syamlal-O’

k -

Brien

0.46

2D-Grid

2D-Grid

0.38

2D-Grid

0.38

2D-Grid 1

0.38

ED

PT

3D-Grid

1 13

0.9

k -

0.46

1 12

K

Syamlal-O’

1 11

ess

RANS

3D-Grid

AC

10

viscosity

Brien

1 9

coefficient,

k -

1 8

RANS

RANS

1 7

coefficient,

CE

6

Frictional

US

3D-Grid 1

3

Specularity

Brien

AN

2

Syamlal-O’

M

1

Restitution

T

Inlet

IP

Grid

CR

Case

Spatial Discretization-Turbulence

l

l 0.9

0.1

None

0.9

0.1

None

ACCEPTED MANUSCRIPT

0.38

1 16

2D-Grid

0.38

1 17

2D-Grid

0.38

1 18

2D-Grid

0.38

1 19

2D-Grid

0.38

1 20

2D-Grid

0.38

1 21

2D-Grid

0.38

1 22

2D-Grid

0.38

0.38

2D-Grid

2D-Grid 2

27

2D-Grid

2D-Grid

2D-Grid 1

k -

Brien

RANS

Syamlal-O’

k -

Brien

RANS

Syamlal-O’

k -

Brien

RANS

Syamlal-O’

k -

Brien

RANS

Syamlal-O’

k -

Brien

RANS

Syamlal-O’

k -

Brien

Syamlal-O’

0.9

0.9

None

0.99

0.1

None

0.95

0.1

None

0.9

0.9

0.99

0.99

0.01

None

0.1

Johnson-et-a l

0.01

Johnson-et-a l

0.9

Johnson-et-a l

0.9

1

Johnson-et-a

k -

Brien

RANS

Gidaspow

0.9

0.1

None

Syamlal-O’

0.9

0.1

None

0.9

0.1

None

0.9

1

Johnson-et-a

Laminar

l

Brien Syamlal-O’

k -

Brien

0.38

RANS

Syamlal-O’

k -

Brien

RANS

Syamlal-O’

k -

Brien

RANS

Syamlal-O’

k -

Brien

RANS

Syamlal-O’

k -

Brien

0.38

0.38

2 29

Syamlal-O’

RANS

0.38

2 28

RANS

PT

0.38

AC

2 26

Brien

CE

2D-Grid 2

25

k -

None

k -

2 24

Syamlal-O’

ED

2D-Grid

RANS

RANS

1 23

Brien

0.5

T

2D-Grid

k -

0.9

IP

15

Syamlal-O’

CR

1

RANS

US

0.38

AN

2D-Grid

M

14

0.46

l 0.9

0.01

Johnson-et-a l

0.9

0.1

Johnson-et-a l

0.9

0.1

Johnson-et-a l

ACCEPTED MANUSCRIPT 30

2D-Grid

0.46

1 31

2D-Grid

0.46

2

RANS

Syamlal-O’

k -

Brien

RANS

Syamlal-O’

k -

Brien

0.9

0.5

None

0.9

0.5

None

FIGURES

IP

7

T

Table 4: Parameters of simulation cases performed in this study.

CR

List of Figures

Figure 1: Graphical explanation of (a) – the restitution coefficient ess and (b) – the specularity

US

coefficient K .

AN

Figure 2: Schematic of a fluidized bed reactor, adapted from [14]

M

Figure 3: Simulation errors for 3D and 2D models.

ED

Figure 4: 3D-RANS simulations. Snapshots at different time instances of (a) – solid volume fraction and (b) – turbulent viscosity ratio. Simulation parameters: drag function is modelled by

PT

Syamlal–O’Brien sub-model; no frictional viscosity sub-model; specularity coefficient K = 0.1

CE

; restitution coefficient ess = 0.9 ; inlet velocity U = 0.38 m/s.

Figure 5: 3D-RANS simulations. Time-averaged solid volume fraction. Simulation parameters:

AC

drag function is modelled by Syamlal–O’Brien sub-model; no frictional viscosity sub-model; specularity coefficient K = 0.1 ; restitution coefficient ess = 0.9 ; inlet velocity U = 0.38 m/s.

Figure 6: Time history of the volume-averaged solid velocity ( U s ) predicted numerically using 2D RANS and 3D RANS simulations with the restitution coefficient of ess = 0.9 , the drag function is modelled by Syamlal–O’Brien sub-model; inlet velocity U = 0.38 m/s

Figure 7: 3D-RANS simulations. Experimental [14] and numerical profiles of time-averaged void

ACCEPTED MANUSCRIPT (gas phase) along x  axis at z = 0.2 m. Simulation parameters: the drag function is modelled by Syamlal–O’Brien sub-model.(a) – U = 0.38 m/s; (b) – U = 0.46 m/s; (c) - comparison between RANS and laminar models for U = 0.38 m/s and K = 0.1 , ess = 0.9 . Figure 8: 2D simulations on Grid1, 56  300 . Solid volume fraction at different time instances

T

obtained as a results of unsteady (a) – Laminar and (b) – RANS simulations. Simulation

IP

parameters: drag function is modelled by Syamlal-O’Brien sub-model; no frictional viscosity

CR

sub-model; specularity coefficient K = 0.1 ; restitution coefficient ess = 0.9 ; inlet velocity

U = 0.38 m/s.

US

Figure 9: 2D simulations on Grid2, 112  600 . Solid volume fraction at different time instances obtained as a results of unsteady (a) – Laminar and (b) – RANS simulations. Simulation

AN

parameters: drag function is modelled by Syamlal-O’Brien sub-model; no frictional viscosity sub-model; no frictional viscosity sub-model; specularity coefficient K = 0.1 ; restitution

M

coefficient ess = 0.9 ; inlet velocity U = 0.38 m/s.

ED

Figure 10: 2D simulation on Grid1, 56  300 . Time-averaged solid volume fraction calculated using unsteady (a1) – Laminar and (a2) – RANS E-E model. Time-averaged gas velocity vector

PT

calculated using unsteady (b1) – Laminar and (b2) – RANS E-E model. Simulation parameters:

CE

drag function is modelled by Syamlal-O’Brien sub-model; no frictional viscosity sub-model; specularity coefficient K = 0.1 ; restitution coefficient ess = 0.9 ; inlet velocity U = 0.38 m/s.

AC

Figure 11: 2D simulation on Grid2, 112  600 . Time-averaged solid volume fraction calculated using unsteady (a1) – Laminar and (a2) – RANS E-E model. Time-averaged gas velocity vector calculated using unsteady (b1) – Laminar and (b2) – RANS E-E model. Case parameters: drag function is modelled by Syamlal-O’Brien sub-model; no frictional viscosity sub-model; specularity coefficient K = 0.1 ; restitution coefficient ess = 0.9 ; inlet velocity U = 0.38 m/s. Figure 12: 2D-RANS simulations. Turbulent viscosity ratio ( t / 0 ) at different time instances. Simulation parameters: drag function is modelled using Syamlal O’Brien sub- model; no frictional

ACCEPTED MANUSCRIPT viscosity sub-model; specularity coefficient: K = 0.1 ; restitution coefficient ess = 0.9 ; inlet velocity U = 0.38 m/s. (a) – Grid1 and (b) – Grid2. Figure 13: Time-averaged gas velocity ( U g ) as a function of height at x = 0.14 m (a)  2D model with frictional viscosity and (b) – 3D model. Simulation parameters: the drag function is

T

modelled by Syamlal O’Brien sub-model; the frictional visocisty is modelled by Johnson-et-al

IP

sub-model; specularity coefficient K = 0.1 ; restitution coefficient ess = 0.9 ; inlet velocity

CR

U = 0.38 m/s.

US

Figure 14: Volume-averaged solid velocity ( U s ) as a function of time for (a) – 2D model with frictional viscosity and (b) – 3D model. Simulation parameters: the drag function is modelled by

AN

Syamlal O’Brien sub-model; the frictional visocisty is modelled by Johnson-et-al sub-model;

M

specularity coefficient K = 0.1 ; restitution coefficient ess = 0.9 ; inlet velocity U = 0.38 m/s.

Figure 15: 2D-RANS simulations. Experimental [14] and numerical profiles of time-averaged

ED

void (gas phase) along x  axis at z = 0.2 m for different values of (a) – restitution coefficient

ess and (b) – specularity coefficient K . Simulation parameters: the drag function is modelled by

PT

Syamlal O’Brien sub-model; no frictional viscosity sub-model; specularity coefficient K = 0.1 ;

CE

restitution coefficient ess = 0.9 ; inlet velocity U = 0.38 m/s.

Figure 16: Time history of the volume-averaged solid velocity ( U s ) predicted numerically using

AC

3D RANS simulations. Simulation parameters: the drag function is modelled by Syamlal O’Brien sub-model; no frictional viscosity sub-model; restitution coefficient ess = 0.9 ; inlet velocity

U = 0.38 m/s. Figure 17: 2D simulations. Experimental [14] and numerical profiles of time-averaged void fraction at z = 0.2 m. Simulation parameters: specularity coefficient K = 0.1 ; restitution coefficient ess = 0.9 ; inlet velocity U = 0.38 m/s (a) – Grid1, and (b) – Grid2.

ACCEPTED MANUSCRIPT

AC

CE

PT

ED

M

AN

US

CR

IP

T

Highlights  Simulations of fluidized bed were performed in ANSYS FLUENT using Euler-Euler model  The results were validated against reference experimental data  3D simulations give more physically realistic results as opposed to 2D simulations  Turbulent flow model should be used to simulate gasification and combustion  The choice of simulation parameters significantly affects numerical predictions

Figure 1

Figure 2

Figure 3

Figure 4

Figure 5

Figure 6

Figure 7

Figure 8

Figure 9

Figure 10

Figure 11

Figure 12

Figure 13

Figure 14

Figure 15

Figure 16

Figure 17