Document [original]

Computational Particle Mechanics

https://doi.org/10.1007/s40571-023-00671-1

Comparison of sub-grid drag laws for modeling fluidized beds with the

coarse grain DEM–CFD approach

Janna Grabowski1·Nico Jurtz1·Viktor Brandt2·Harald Kruggel-Emden2·Matthias Kraume1

Received: 7 July 2023 / Revised: 1 September 2023 / Accepted: 22 September 2023

Abstract

Fluidized particulate systems can be well described by coupling the discrete element method (DEM) with computational fluid

dynamics (CFD). However, the simulations are computationally very demanding. The computational demand is drastically

reduced by applying the coarse grain (CG) approach, where several particles are summarized into larger grains. Scaling rules

are applied to the dominant forces to obtain precise solutions. However, with growing grain size, an adequate representation

of the interaction forces and, thus, representation of sub-grid effects such as bubble and cluster formation in the fluidized

particulate system becomes challenging. As a result, particle drag can be overestimated, leading to an increase in average

particle height. In this work, limitations of the system-to-grain ratio are identified but also a dependency on system width. To

address this issue, sub-grid drag models are often applied to increase the accuracy of simulations. Nonetheless, the sub-grid

models tend to have an ad hoc fitting, and thorough testing of the system configurations is often missing. Here, five different

sub-grid drag models are compared and tested on fluidized bed systems with different Geldart group particles, fluidization

velocity, and system-to-grain diameter ratios.

Keywords Fluidized bed ·Computational fluid dynamics ·Discrete element method ·Sub-grid drag model ·Coarse grain

particle method

List of symbols

Abbreviations

CFB Continuous fluidized bed

CFD Computational fluid dynamics

CG Coarse grain

DEM Discrete element method

DNS Direct numerical simulation

DPVM Divided particle volume method

EMMS Energy minimization multi-scale

FFT Fast Fourier transformation

norm. Normal

PSD Power spectral density

s-s Steady-state

tan. Tangential

TFM Two fluid model

BJanna Grabowski

janna.grabo[email protected]

1Chair of Chemical and Process Engineering, Technische

Universität Berlin, Ackerstr 76, 13355 Berlin, Germany

2Chair of Mechanical Process Engineering and Solids

Processing, Technische Universität Berlin, Ernst-Reuter-Platz

1, 10587 Berlin, Germany

Greek letters

βInterphase momentum transfer coefficient, kg/(m2s)

δOverlap, m

γDamping coefficient, Ns/m

μFriction coefficient

μfFluid viscosity, Pa s

ωAngular velocity, s−1

ρDensity, kg/m3

εVoid fraction

Latin symbols

UdSuperficial gas velocity of dense phase, m/s

FcContact force, N

FdDrag force, N

FhTotal hydrodynamic force, N

F∇pPressure gradient force, N

gAcceleration due to gravity, m/s2

TTorque, Nm

ufFluid velocity, m/s

uslip Slip velocity, m/s

vpParticle velocity, m/s

fFilter size, m

123

Computational Particle Mechanics

f∗Dimensionless filter size

CdDrag coefficient

DSystem diameter, m

dParticle or grain diameter, m

dx CFD cell size, m

eRestitution coefficient

hkAxial position of center of cell k,m

IpMoment of inertia, m2kg

kSpring constant, N/m

lCoarse grain factor

LcCharacteristic cluster length scale, m

Ncells Number of cells

ncg Coarse grain number

uin Fluid inlet velocity, m/s

umf Minimum fluidization velocity, m/s

u∗

slip Filtered, dimensionless slip velocity

VkVolume of cell k, m3

vtTerminal particle velocity, m/s

∇pPressure gradient, Pa/m

Ar Archimedes number

hBed height, m

He Heterogeneity index

mMass, kg

NRatio of system diameter to particle/grain size

Re Reynolds number

VVolume, m3

Subscripts

ave Average

f Fluid

gGrain

het Heterogeneous

init Initial

nNormal

p Particle

ref Reference

t Tangential

1 Introduction

Fluidized beds play an essential role in the chemical and

process industry [1]. An upward-directed gas stream flu-

idizes a bulk of particles in a fluidized bed. Their good heat

and mass transfer characteristics make it an efficient sys-

tem for multiple chemical or physical applications [2,3].

However, such systems are inherently unstable, and effects

occur on multiple time and length scales. On a microscopic

level, particle–particle and particle–wall interactions, as well

as the interaction of particles with the fluid, are dominant.

On the mesoscale level, bubble- and cluster formation domi-

nate due to local fluid flow effects. These strongly affect the

macroscale, e.g., the fluidization regime and mixing charac-

teristics [4].

Many approaches have been proposed to model particle–

fluid multiphase systems to account for the different time

and length scale phenomena [5].

Direct numerical simulations (DNS) can be applied on a

microscopic scale. Here, fluid dynamics is solved by compu-

tational fluid dynamics methods such as Lattice Boltzmann

(LBM) or finite volume method (FVM). The particle sur-

face can be resolved using the immersed-boundary method

[6–8]. Particle movement is solved with the discrete element

method. The methods are computationally very demanding

and limited to small-scale systems [9].

Mesoscale effects can be resolved using either the two-fluid-

model (TFM) or coupling DEM and particle-unresolved CFD

(DEM–CFD). TFM is an Euler–Euler approach in which

particle–particle, particle–wall, and particle–fluid interac-

tions are modeled by closure conditions such as the kinetic

theory of granular flow [10]. Much research is focused

on developing accurate closure models, as shown by [11–

13]. DEM–CFD is a Lagrange–Euler method where particle

movement is tracked by solving Newton’s second law.

Particle–particle and particle–wall interactions are described

using contact models such as the Hertzian model or the linear-

spring-dashpot model [14]. Closures are only necessary to

model the fluid-particle interaction with appropriate drag

laws. However, the well-resolved TFM needs a grid resolu-

tion with a cell size less than 10 particle diameters [15], which

is not feasible for industrial applications. Also, the DEM–

CFD approach is computationally very demanding and not

feasible for industrial systems with particle numbers ranging

from 108to 1012 [5].

A common way to decrease the computational demand is to

reduce the cell number by using a coarse mesh. This can be

done for both TFM and DEM–CFD. In the TFM, both the

description of the solid and fluid phases is affected. In the

DEM–CFD only the fluid phase is implicated by a mesh

coarsening. In the DEM, the computational demand with

regard to the solid phase can be decreased by reducing the

number of particles following the coarse grain (CG) approach

[16]. Here, several original particles are lumped into a larger

grain under the assumption that the individual movement of

a particle is not crucial to the fluidized bed mesoscale char-

acteristics as shown by [17]. The CG approach gives the

possibility to simulate large-scale reactors in a realistic time

frame. Some large-scale examples are given in [18–20].

However, if the cells as part of a TFM or CG-CFD-DEM are

too coarse in combination with a too strong coarse graining

for the latter, it can lead to an inadequate representation of

the mesoscale characteristics. A particle moving in a par-

ticle cluster experiences less drag than individual particles

[1]. Neglecting this phenomenon in coarse simulations, thus

assuming the drag in a grain is similar to the sum of the

123

Computational Particle Mechanics

particles in the grain, where each particle exerts drag as a

single particle, can lead to an overestimation of drag force

and, thus, an erroneous estimation of bed expansion [21–

23]. To address the aforementioned issue, much research

has been conducted on so-called sub-grid drag models that

account for the mesoscale effects while still using a coarse

mesh or coarse mesh and coarse grains. Sub-grid drag mod-

els can be classified into empirical, structure-based, energy

minimization-based, filtered, and artificial neural network-

based models [24]. The sub-grid drag model can be related to

empirical correlations, for example, the O’Brien and Syam-

lal drag law, which is based on an air-FCC system with solid

circulation fluxes [25]. Wang et al. [26] developed a model

based on the assumption that for Geldart B and D particles,

interphase drag, e.g., drag between the solid and fluid phase,

only occurs in a particle-dense phase. In contrast, the dilute

phase, e.g., bubbles, remains free of particles and free of inter-

phase drag. This leads to a two-phase structure of the system

for the drag estimation, which belongs to the category of a

structure-based drag law.

A large and well-investigated drag law group is based on

the concept of energy minimization of particle transport and

suspension [27,28]. The energy minimization multi-scale

(EMMS) model has been expanded to account for bubbly

and circulating fluidized beds recently [29,30].

Lastly, filtered drag correlations are often applied. Here, a

correction factor is introduced, obtained by regressing results

of finely resolved TFM or DEM–CFD simulations [31–33].

The filtered drag models can tend to ad hoc regressed results,

so applicability to other system configurations always needs

evaluation. A recent approach in the development of filtered

drag laws involves leveraging artificial neural networks for

regression purposes [34,35]. Lu et al. [35] performed a com-

parative analysis between filtered drag laws derived with

artificial neural network (ANN) and the nonlinear regres-

sion technique. In this case, the filtered drag law generated

through nonlinear regression showed more accurate results.

Moreover, attempts to couple the ANN method to the EMMS

model were successfully made to extend EMMS to multiple

macroscale markers such as operating conditions, e.g., gas

inlet velocity, system width, or even material properties that

are averaged over the cross-sectional area of the fluidized bed

system [36,37].

The heterogeneous drag laws are, in most cases, devel-

oped for the application of coarse TFM simulations. Few

researches have been conducted on the application of het-

erogeneous drag laws in the coarse grain DEM approach

framework. Lu et al. [38] demonstrated the general applica-

bility of the EMMS model with coarse grains and showed

good agreement with experimental data. Moreover, Jurtz et

al. [39] investigated the EMMS approach for a bubbling flu-

idized bed with Geldart class A particles, which improved

the results compared to a classical homogeneous model. The

EMMS model has also been applied to DEM simulations [40]

and for an industrial example of a polymerization process

using coarse grains [41], which was also in good agreement

with experimental data. Radl et al. [31] extended their filtered

drag law to the CG approach by introducing a factor reducing

the drag on the grains with increasing grain size. However,

they stated that their model might be ad hoc regressed, and

it is only tested for coarse grain factors up to 8.

As can be seen, research about the applicability of such mod-

els to the CG-DEM–CFD approach has been demonstrated

for EMMS and filtered drag models. However, a detailed

analysis of the applicability of sub-grid drag models to CG-

DEM–CFD simulations is still missing. This work compares

popular drag models using the CG-DEM–CFD approach for

fluidized beds. Often the drag models are limited to specific

operating conditions of the system. Therefore, limitations of

the drag laws concerning fluidization regime, system width,

and particle types (Geldart class B and D) are analyzed while

changing the coarse grain factor, thus, the number of grains

in the system. Geldart class A particles play an essential role

in fluidized bed engineering and will be the focus of future

work.

2 Numerical methodology

The DEM is used to solve the movement of particles. Within

this method, the soft sphere approach is employed with the

linear spring-dashpot model to model the interaction between

particles. For fluid movement, computational fluid dynamics

(CFD) is used based on fluid phase continuity and momentum

conservation. The DEM and CFD are coupled and solved

as part of the software MFiX [42,43]. As MFiX does not

provide for the sub-grid drag laws considered here, MFiX

is extended using the usr_drag module provided by MFiX.

Regarding the DEM–CFD, a full description of the model

and applied equations is given in Sect. A.1.

2.1 Coarse grain DEM

Coarse graining is performed by summarizing several orig-

inal particles into larger grains. The coarse grain number

ncg =np/ngis defined as the number of particles lumped

into a grain. The coarse grain factor, also known as grain-to-

particle size ratio, is defined as:

l=dg/dp,(1)

where dgis the grain diameter and dpis the particle diameter.

Contact forces in CG-DEM are calculated by assuming com-

pact grains, so the total volume and the total mass are

constant. A rescaling of the interaction forces is needed

to represent the original particles properly. Several scaling

123

Computational Particle Mechanics

approaches have been developed recently. For our simulation

system, the scaling method by [44] combined with [45]was

proven one of the most accurate looking at multiple output

parameters in our previous work [46] and it is applied here.

The stiffness components in normal and tangential directions

are scaled with the coarse grain number:

kn,g=kn,p·ncg,(2)

kt,g=kt,p·ncg,(3)

whereas, for the damping components, the ratio of normal to

tangential damping is kept constant. The restitution coeffi-

cient is modified to represent intra-grain collisions such as:

eg=e√ncg

p.(4)

Lastly, when sliding occurs, dissipation is conservated by

scaling the friction coefficient μas follows

μg=μp

√l.(5)

As hydrodynamic forces, pressure gradient force F∇p

pand

drag force Fd

pare considered. The pressure gradient force

acting on the grain is assumed to be similar to the sum of

particles inside the grain:

F∇p

g=ncgF∇p

p=Vg∇p(6)

with Vgas the volume of the grain. In our previous work [46],

it was found that the scaling of drag force in a similar way

gives the most accurate results, such as:

g=ncgFd

p=Vg

1−εf

(uf−vg), (7)

where εfis the void fraction, ufand vgare the fluid and grain

velocity, respectively, and βis the interphase momentum

transfer coefficient. Drag models depending on slip velocities

and void fractions are applied to calculate this coefficient. In

this study, the drag law developed by Gidaspow [47] is used,

which belongs to the class of homogeneous drag models. It

is a combination of the Ergun drag law [48] for dense regions

and Wen-Yu drag law [49] for dilute regions:

β=⎧

⎨

⎩

150 (1−εf)2

εf

μf

p+1.75(1−εf)ρ f|(uf−vp)|

dp,εf<0.8

4dpCdεf(1−εf)ρ f|(uf−vp)|ε−2.65

f,εf>0.8,

(8)

where μfis the fluid viscosity, ρfis the fluid density and vp

is the particle velocity. The drag coefficient Cdis calculated

as:

Cd=24

Rep(1+0.15Rep)0.687,Rep<=1000

0.44,Rep>1000,(9)

where Repis the particle Reynolds number, which is defined

as:

Rep=|uf−vp|εfdp

μf

.(10)

However, for coarse simulations, either by applying a larger

coarse grain factor or by applying a larger grid size, the accu-

rate representation of sub-grid behavior becomes challeng-

ing. Therefore, modification of the interphase momentum

transfer coefficient is performed by introducing a heterogene-

ity index that accounts for the mesoscale effects:

He =βHet

βWen-Yu

.(11)

Here, βWen-Yu is the interphase momentum transfer coeffi-

cient obtained by the Wen-Yu drag law. βHet is the interphase

momentum transfer coefficient obtained by heterogeneous

sub-grid drag laws. The sub-grid drag laws will be explained

in the following section.

2.2 Sub-grid drag laws

This subsection elaborates on the sub-grid drag laws used

in this study. All necessary equations used to implement the

drag laws in the software MFiX are given in Sect. A.2.The

structure-based drag law is explained first, followed by the

EMMS drag law group and the two filtered drag laws.

Structure-based drag law: As a structure-based drag law,

the drag law by Wang et al. [26] is applied. The modi-

fication by Wang et al. to the homogeneous drag model

is based on the assumption that the drag between solid

and fluid only occurs in a particle-dense phase. For Gel-

dart class B and D beds, bubbles remain primarily free

of particles, which justifies the assumption for the drag

for these cases. In the model, the drag force in each

CFD cell is adjusted to only account for the fraction

of the cell occupied by the particle-dense phase. The

fraction of the particle-dense phase is calculated based

on a Reynolds–Archimedes correlation. The interested

reader is referred to Wang et al. [26] for detailed infor-

mation. The model was initially designed for the TFM.

The modified model improved results compared to the

unmodified drag law looking at a fluidized bed system

at the pilot and at the industrial scale [26]. For the CG-

DEM–CFD approach, the model has not been tested yet.

Nevertheless, the underlying assumptions are also valid

123

Computational Particle Mechanics

for a Lagrange–Euler framework, which is why this work

investigates the model.

EMMS group: Two EMMS drag laws, the EMMS/Matrix

[27,28], and the EMMS/steady-state model [29,30],

are investigated in this work. The EMMS/Matrix drag

law, initially developed for flow in risers, accounts for

the heterogeneous behavior by establishing transport

equations of the particle flow in the dense, particle-rich

phase, the dilute, fluid-rich phase, and the interaction

of dense and dilute phase of the fluidized bed. This set

of equations leads to an unclosed system of equations.

Therefore, a stability criterion is introduced to minimize

the energy consumption for suspending and transport-

ing the particles in the flow [27,28]. The solution to the

minimization problem is an interphase momentum trans-

fer coefficient for a specific range of slip velocities and

void fractions. From this, the heterogeneity index can be

calculated with Eq.11. In this work, the EMMS/Matrix

minimization problem is solved before the DEM–CFD

simulation. From the model, a table is generated with

the heterogeneity index for a specific range of void frac-

tions (0.42 <ε

f<0.994)and slip velocities (for

Geldart D: 0.001 m/s <uslip <1.4 m/s, for Geldart B:

0.001 m/s <uslip <2.77 m/s). During the simulation,

the heterogeneity index in each CFD cell is obtained by

bilinear interpolation of the void fraction and slip veloc-

ity in the CFD cell to the values in the table.

The EMMS model was extended to lower throughput

beds such as circulating fluidized beds and bubbling

fluidized beds by solving a steady-state EMMS model

[29,30]. Here, the heterogeneity index depends only

on the according void fraction in the CFD cell. The

EMMS/steady-state model is incorporated into the soft-

ware (MFiX) using a correlation, given in Sect. A.2.

Filtered drag laws: Radl et al. [31] set up a filtered drag

model based on well-resolved Lagrange–Euler simula-

tions. The well-resolved Lagrange–Euler simulations of a

3D periodic domain were performed using a constant par-

ticle diameter dpof 75 µm and varied average porosity.

With these simulations, averaged interphase momentum

transfer coefficients were determined for various filter

sizes franging from 3.3dpto 30dp, which can be used

for coarse simulations. Moreover, a new characteristic

cluster length scale Lcis established for coarse simula-

tions with different particle types based on the Froude

number. The filtered drag law encompasses both fluid

coarsening, which refers to the coarsening of the CFD

cell size, and coarse graining. A correction factor is intro-

duced for coarse grain DEM to decrease the interphase

momentum transfer coefficient for increasing grain size.

However, their model is only tested for l<8, so inaccu-

racies may occur for larger grain diameters and different

Fig. 1 Heterogeneity index trends with changing void fraction at a

constant slip velocity of uslip =0.98 m/s for the applied drag laws

system configurations.

Similar to Radl et al., the filtered drag model by Lu

et al. [35] is based on averaging well-resolved DEM–

CFD simulations of a fluidized bed. It is a revised

version of the filtered drag law by Sarkar et al. [50]. The

original drag law was developed for TFM simulations.

Now, the revision allows using the drag law for DEM–

CFD simulations, which makes it more applicable for

CG-DEM–CFD. Moreover, by introducing a macroscale

marker for the gas inlet velocity, the model is applica-

ble to systems with a wider range of gas inlet velocities.

It was shown that coarse grain simulations performed

with filtered drag laws derived from DEM–CFD simula-

tions yield more accurate results than filtered drag models

derived from TFM simulations [35,51]. Similarly to the

Radl drag law, the heterogeneity index is calculated with

regression parameters that differ with void fraction and

slip velocity. More information is available in [35].

A comparison of the heterogeneity index as a function of

the void fraction for the different sub-grid drag laws is given

in Fig.1for a constant slip velocity of 0.98 m/s for Geldart

B particles and a gas inlet velocity of vin =1.7 m/s. For the

regression based models, the ratio of the filter size to charac-

teristic cluster length scale f/Lcequals 1.08. The Wang

drag law [26] decreases drag compared to the Gidaspow drag

law until a void fraction of εf=0.55, leading to a hetero-

geneity index below 1. From there, the heterogeneity index

reaches values up to 4. The EMMS/Matrix model, as origi-

nally established for fast risers, lowers the drag for the void

fraction area 0.42 <ε

f<0.99 between 1 and 2 orders

of magnitude. No correction is performed for void fractions

below 0.42 and above 0.99. The EMMS/steady-state model

reduces the drag for 0.45 <ε

f<0.85 up to one order of

magnitude, but a drag correction is not conducted for dense

and dilute areas. The drag law by Radl et al. shows a reduction

of He =1toHe=0.7 up to a void fraction εf=0.6. Subse-

quently, the heterogeneity index slightly increases again but

123

Computational Particle Mechanics

Table 1 Simulation parameters

Unit Geldart D simulations Geldart B simulations

Particle diameter mm 1 0.5

Particle density kg/m32500 1500

Spring constant (normal) N/m 2 ×1041×103

Spring tan./norm. ratio – 0.32 0.32

Damping tan./norm. ratio – 0.1 0.1

Friction coefficient – 0.1 0.1

Restitution coefficient – 0.97 0.97

Fluid density kg/m31.2 1.2

Fluid viscosity Pa s 1.8×10−51.8×10−5

Fluid inlet velocity m/s 1.35/2.75 0.5/1.7

Min. fluidization velocity m/s 0.53 0.12

System height mm 400 200

System diameter mm 100/75/50 50/37.5/25

Initial average particle height mm 25 12.5

Particle number – 4 ×105/2.25 ×105/1 ×1054×105/2.25 ×105/1 ×105

is always below 1. The drag law by Lu et al. decreases drag

significantly for a void fraction of εf=0.42 and εf=0.99,

but only minor drag reduction is conducted for a void frac-

tion range of εf=0.6−0.8. No correction is performed

for void fractions below 0.42 and above 0.99, similarly to

the EMMS/Matrix model. Note that the heterogeneity index

calculated with the drag laws by Radl et al., Lu et al., and the

EMMS/Matrix law also depend on the slip velocity. Thus,

the actual behavior can be different based on the reached slip

velocities in the system.

3 Simulation design

The simulations are designed to assess the drag laws’

application range thoroughly. For this, a cylindrical system

configuration is chosen with three different system widths

to assess whether a limitation on the coarse grain factor

or the ratio of system diameter to particle/grain diameter

N exists. Moreover, two gas inlet velocities are chosen to

account for two fluidization states. Monodisperse particle

beds are generated with Geldart B or Geldart D particles to

assess the effectiveness of the sub-grid drag laws for dif-

ferent particle types. Simulations are performed with the

open-source software MFiX [42,43]. The properties for the

unscaled benchmark particle systems are given in Table 1.

For the particle properties, the properties of the particle–

particle interactions are given. The properties of particle–wall

interactions are the same, except for the friction coefficient,

which is reduced from 0.1 to 0.09. Rolling friction is not con-

sidered because appropriate scaling rules of rolling friction

for the coarse grain approach are still a subject of current

research [16]. The benchmark simulations are then repeated

with coarse grain factors l=2,4,6,8 and 10, applying the

scaling rules explained in Sect. 4.1.

The system width is reduced from D=100 to D=75 mm

and D=50 mm for Geldart D particles and from D=50

to D=37.5 mm and D=25 mm for Geldart B particles.

This leads to a constant system-to-particle ratio Npof 100

to 50 and 25, respectively. A constant initial bed height is

ensured by reducing the particle number from 4 ×105to

2.25 ×105and 1 ×105particles, respectively. The spring

constant is adjusted for Geldart B and Geldart D type particles

to keep the maximum overlap between particles or particles

and wall below 1% as recommended by [52] but ensure a

reasonable DEM time step. The ratio between particle col-

lision time and DEM time step is set to 25. The CFD time

step is set constant to 1 ×10−4s. The CFD cell size is set

to dx =2dpfor unscaled simulations and dx =2dgfor

coarse grain simulations in each direction. The cells are cut

into smaller cells at the boundaries to resolve the cylindrical

geometry of the domain using the cut-cell method. There-

fore, the divided particle volume method (DPVM) is used

for mapping Lagrangian information on the Eulerian mesh

to avoid numerical instabilities. Traditionally, void fraction

and other interphase quantities are assigned to the CFD cell,

in which the centroid of the particle is, which is called the

particle centroid method [53]. For DPVM, these quantities

are also distributed to the neighboring CFD cells based on

a distributing weight function, which allows using smaller

cell-to-particle ratios [43,53]. To map the gas velocity back

to the particle’s position, two steps are performed in MFiX.

At first, the gas velocities are mapped to an interpolation grid.

At second, a second-order accurate Lagrange polynomial is

123

Computational Particle Mechanics

Fig. 2 Geometries for the three different system widths and boundary

conditions

used to interpolate the gas velocity from the interpolation

grid to the particle’s position [54]. The gas flows through

the domain using a velocity boundary condition at the bot-

tom and a pressure outlet at the top, as illustrated in Fig.2.

No-slip boundaries are set for the walls of the fluidized bed.

The physical time simulated is 20 s. The gas inlet velocity is

ramped in the first second of the simulation to ensure smooth

fluidization.

For evaluation of the simulation accuracy, the normalized

average particle or grain height is calculated. Therefore, at

first, the average particle or grain height h<p/g>for each

time frame is calculated using the formula by [55]:

h<p/g>=Ncells

k=0(1−εf,k)hkVk

Ncells

k=0(1−εf,k)Vk

,(12)

where hkis the axial position of the center of cell k and Vk

is the volume of the cell k. At second, the average particle or

grain height is then normalized by dividing the initial average

particle height h<p,init>of the static bed given in Table 1.

The general simulation procedure is as follows: In the first

step, unscaled DEM–CFD simulations with the Gidaspow

drag law are performed for each system configuration, which

is used as a reference. Three different system diameters

are investigated to analyze the limitation of the system-to-

particle/grain diameter ratio N. In the second step, each

configuration is repeatedly simulated using coarse grain fac-

tors between l=2–10 with the homogeneous Gidaspow

drag law. The simulations with a system-to-particle ratio of

Np=50 are only repeated until l=8 to ensure at least

3 CFD cells in each direction and maintain a cell size of

dx =2dp/g. The limitations of resolving sub-grid phenom-

ena as a function of the system-to-grain diameter ratio or

the coarse grain factor are analyzed and determined. In the

third step, the simulations with large coarse grain factors of

l=6−10 and l=6−8forNp=50 are repeated with

the heterogeneous drag laws defined in Sect. 2.2. Most sub-

grid drag models are optimized on a specific particle type

and fluidization regime. Geldart group B and Geldart group

D particles are analyzed to assess applicability for different

particle types. The simulations are conducted in bubbling and

slugging to turbulent regimes to cover two different types of

fluidization.

4 Results and discussion

As described above, this section is divided into two sub-

sections to analyze and discuss the impact of drag laws on

the accuracy of the coarse grain approach. Firstly, the cases

with the homogeneous Gidaspow drag law are analyzed for

various coarse grain factors. Secondly, the cases showing

significant deviation from the unscaled benchmark cases are

then analyzed using the different sub-grid drag laws. The goal

is an improvement of the results for all tested configurations.

4.1 Effect of coarse graining using the Gidaspow

drag law

In Fig.3, the normalized average particle or grain height for

Geldart class B (Fig.3a) and D (Fig.3b) particles is shown

by using the Gidaspow drag law under varying system-to-

particle/grain diameter ratios N. The ratio N is varied by

changing the system diameter and the particle/grain diam-

eter. Here, the normalized average particle or grain height

is temporally averaged over the time frame of t=5−20 s,

which results in one value for each simulation.

Looking at the Geldart class B systems, almost no impact of

coarse graining with respect to the normalized particle height

can be observed until a ratio of N<12.5 for both inlet

velocities, except for the system with a ratio of Np=50,

where the bed height is increasing only at a ratio lower than

N<9.4. A ratio of N=12.5 corresponds to a maximum

coarse grain factor of l=8 in the wide system and l=4

in the narrow system, as seen in the color bar, representing

the coarse grain factor l.BelowN<12.5 and N<9.4,

a relatively sharp increase in particle height is observable.

Interestingly also is the influence of absolute system width.

For similar ratios N, the normalized particle height in the

wide system only increases from 4.04 to 4.27 even for coarse

grain factors of l=10 (N=10), whereas for the width of

Np=50 and Np=75, the bed height rises sharply from

4.04 and 4.17 in the benchmark simulations to 4.73 and 4.93

for l=10 (N=7.5) and l=8(N=6.25).

123

Computational Particle Mechanics

Fig. 3 Normalized particle height for coarse grain factors between 2

and 10 for Geldart B particles (a) and Geldart D particles (b)atthe

two gas inlet velocities for the three system widths. The symbols mark

the different system widths, unfilled at the low gas inlet velocity and

filled at the high gas inlet velocity. The different colors of the symbols

represent l

For Geldart D class systems, at the large inlet velocity of

vin =5.2umf, the normalized particle height is constant for

ratios N until N=25, with one exception for Np=75.

There, at N=37.5, the value is already overpredicted,

meaning the deviation to the benchmark case is more than

5%. Similar, to the Geldart class B, there is also an absolute

influence of the system width as a narrower system width is

leading to larger normalized particle heights at the same ratio

N. For the low gas inlet velocity of vin =2.6umf, the normal-

ized particle height is first underestimated and then increased

for larger grain sizes. Therefore, no limitation regarding the

coarse grain factor or ratio Ncan be concluded for the low

gas inlet velocity for Geldart class D particle systems.

As shown in Fig.3, an overprediction of the particle height

appears already at larger Nfor both particle types when using

a larger gas inlet velocity. The Geldart class D particles are

less elevated because the fluidization regimes differ depend-

ing on particle type. However, the deviation with decreasing

ratio Nis slightly larger for Geldart D particles. For Geldart

class D, the lowest normalized particle heights are reached

at around 2.6−2.8 for benchmark cases, but reach values up

to 4 for low ratios N. For Geldart class B, the increase is

from 4 for benchmark cases to 5 for low ratios N. The rea-

son for the larger error using Geldart class D particles might

be explainable by two effects.

Firstly, at the large inlet velocity, different fluidization

regimes are reached for Geldart B particles and Geldart D

particles. A good indicator for the fluidization regime is the

intensity of pressure fluctuations in the fluidized bed [1,56].

Therefore, a fast Fourier transformation (FFT) on the pres-

sure drop signal is performed. The pressure drop is calculated

between the inlet and the outlet of the fluidized bed system.

A sampling rate of 100 Hz considering data taken from the

time frame of t=5−20 s is chosen. In Fig. 4, the results of

the FFT are shown in the form of the power spectral density

(PSD) plotted against frequency for the unscaled benchmark

cases at the two velocities.

Exemplary, the results from the widest system (Np=100)

are shown. For the simulations with the Geldart class B

particles, the power spectral density always stays below

1500 Pa2/Hz. The main peaks at the low gas inlet velocity are

between 5 and 15 Hz, and the PSD goes up to 250 Pa2/Hz.

For the gas inlet velocity at 14umf, the dominant peak is at

5 Hz with a much higher PSD but narrower width, which is

typical for the slugging regime [1,56]. Looking at the simu-

lations with Geldart D particles, the power spectral density is

one magnitude higher. Interestingly, the simulation at lower

gas inlet velocity shows much higher power spectral den-

sities at specific frequencies than at high gas inlet velocity.

A decrease in the PSD is typically seen when going from a

bubbling or slugging regime to a turbulent regime [1,56].

Therefore, the frequency analysis indicates that the system

with Geldart D particles at the high gas velocity is reaching

the turbulent regime. In contrast, the system with Geldart

B particles is reaching the slugging regime. For the turbu-

lent regime, a higher heterogeneity of the fluidized bed is

expected. This might lead to overprediction in particle height

as the coarse systems cannot resolve the heterogeneous char-

acteristics of the fluidized bed.

Secondly, the difference in absolute gas inlet velocity can

lead to the different behavior of Geldart B to Geldart D grains.

The fluid velocity to achieve the higher fluidization regime

is uin =2.75 m/s for Geldart class D particles, and uin =

1.7 m/s for Geldart class B. This could result in higher slip

velocities being attained by Geldart D particles, potentially

123

Computational Particle Mechanics

Fig. 4 Power spectral density

plots for the unscaled

benchmark cases in the wide

system (Np=100) for Geldart

class B (a)andD(b) particles.

The blue line shows the result at

the low gas inlet velocity, and

the black line at the high gas

inlet velocity

(a) Geldart B (b) Geldart D

amplifying the error in slip velocity calculations due to the

larger cell size.

Moreover, it is also observable that the coarse grain results in

the narrow system Np=50 have a higher deviation from the

benchmark than in the wider systems. Therefore, there might

also be other effects causing the deviation in the narrow sys-

tems for both Geldart types. There can be four reasons why

the effect cannot be corrected. Firstly, the D/d-ratio Ncan

be too narrow to distribute the particles in the radial direction

leading to larger distribution in the axial direction. Secondly,

the portion of particle–wall interaction compared to particle–

particle interaction increases significantly and might lead to

a deterioration of the results. Thirdly, the large grains distort

the typical fluid flow pattern in narrow configurations. Lastly,

due to the inherent mesh coarsening, the adequate resolu-

tion of fluid flow and interaction quantities for drag, such as

porosity and slip velocities, become challenging. For exam-

ple, for l=8 and Np=50, there are only 9 grid points in

the radial direction. In addition, the coarse grid might lead to

erroneous interpolation of the velocities between the no-slip

boundary condition at the cylinder wall to the fluid velocity

at the position of the particle depending on the interpolation

scheme applied. A second-order accurate Lagrangian poly-

nomial is already used. Switching to a higher order, in this

case a fourth-order accurate Lagrangian polynomial did not

enhance the results and lead to an even higher overestimation

of the normalized particle height. Interpolation techniques

exist to reduce the CFD cell-to-particle ratio and might lead

to improved results. However, interpolation errors can appear

and must be assessed for such methods, so they were not

applied in these cases.

Figure5shows snapshots of the wide simulation system

(Np=100) taken at 10 s to give a visual example of

the change of the fluidized bed dynamics while applying

coarse graining. Exemplary, the simulation results for par-

ticle velocity for Geldart B particles with an inlet velocity of

uin =14.2umf are shown. For the benchmark case, l=1,

we can see most of the particles concentrating in the lower

part of the system. Some particles move into the upper part,

indicating a dense region and a dilute region. The separation

Fig. 5 Snapshots taken at 10 s of a Geldart B type simulation with an

inlet velocity of 14umf in the wide system of Np=100. ashows the

benchmark case (l=1),bfor l=4, cfor l=10. On the right-hand

side, the time-averaged void fraction profiles for these cases over the

system height are shown

between the dense and the dilute phases is still observable

for a coarse grain factor of l=4, which reduces the par-

ticle number from 400,000 to 6250. Although particles are

not elevated as high as in the unscaled benchmark case, the

bed average void fraction profile is shifted to a lower void

fraction with increased bed height. This effect is even inten-

sified for the coarse grain factor of l=10, as also already

shown in Fig.3, the normalized particle height is increased.

The snapshot shows that particles are elevated, but differen-

tiation between the dense and dilute regions is not possible

anymore, as one grain holds 1000 particles.

A similar effect is also observable by looking at the contour

plots of the void fraction in Fig.6. Here, exemplary shown for

Geldart class D particles in the wide system (Np=100) at a

gas inlet velocity of uin =5.2umf.Forl=1 and l=2, sev-

eral bubbles and small particle clusters are visible. For l=4

and l=6, the bubbles turn into larger voids. A differentia-

tion between dense and dilute regions becomes challenging

for large coarse grain factors of l=8 and l=10.

To sum up the subsection, the simulations performed with

the Gidaspow drag model show no overprediction of particle

height for low coarse grain factors and thus, larger D/d-ratios

123

Computational Particle Mechanics

Fig. 6 Representative contour plots for void fraction profiles taken at 5

s of the simulation for Geldart D type particles at fast fluidization uin =

5.2umf and a system width of Np=100. The black line represents a

void fraction of ε=0.9. The decreasing resolution of the simulation

leads to erroneous bubble determination and a loss of the ability to

determine particle clusters, which is observable looking at the slices

going from l=4tol=6

N. When using larger coarse grain factors, so low ratios N,

the bed height is overpredicted, probably due to an overesti-

mation of the drag force. This effect is intensified when using

higher gas inlet velocities. At low gas inlet velocities, most

particles concentrate at the bottom of the system and the influ-

ence of drag force on the particles is rather small, which leads

to accurate results or even underprediction of particle height.

However, for larger gas inlet velocities, the heterogeneous

behavior of the fluidized bed can be intensified, which might

be less resolved with the coarse simulations. As the devia-

tion is dependent on the gas inlet velocity, the overestimation

might increase with large gas velocities, looking, for exam-

ple, at circulating fluidized beds or risers. When comparing

the fluidized bed simulations containing Geldart D particles

with the simulations containing Geldart B particles, the ten-

dencies are similar, especially for the low inlet velocities. For

the large inlet velocities, the deviation for Geldart D particles

is more exaggerated than with the Geldart B particles.

4.2 Comparison of sub-grid drag models

Heterogeneous drag models tend to be ad hoc modified,

meaning they can be restricted to the simulation system they

were generated with. Hence, the following requirements are

established. The modified drag model should decrease the

observed deviations of normalized particle height for coarse

grain simulations. The model should work for different parti-

cle types (Geldart class B and D particles are tested), different

fluidization regimes (two inlet velocities are investigated),

and different system-to-grain diameter ratios (three different

system diameters are analyzed).

Fig. 7 Boxplots of the normalized particle height for Geldart B and

Geldart D particles. Here, exemplary shown at the high gas inlet velocity

(for Geldart class D: uin =5.2umf and for class B: uin =14.2umf)in

thewidesystem(Np=100) for a coarse grain factor of l=10

4.2.1 Effect of particle type

Figure 7shows the normalized particle height characteris-

tics for the heterogeneous drag laws compared to the DEM

benchmark for Geldart class D and B particles using the large

inlet velocity and a coarse grain factor of l=10 for the wide

system of Np=100. The unscaled benchmark results are

shown on the left for comparison. The distribution of average

normalized particle height over time (t=5−20 s) is summa-

rized in a boxplot, hence, showing the minimum, maximum,

median, and quartiles of the normalized particle height in a

concise way. Outliers are shown as red symbols. An outlier

is a single data point that is larger than 1.5 times the value of

the upper quartile (Q3) or 1.5 times lower than the value of

the lower quartile (Q1). Therefore, the boxplots not only give

information about average results but also about the dynam-

ics of the fluidized bed system.

For the benchmark, the distribution of particle height shows

a symmetrical trend with a similar distribution of values

below and above the median. For Geldart D, more outliers are

observed as for the Geldart B simulation. These tendencies

are also observed using the Gidaspow model with l=10 but

with a general trend of overestimating normalized particle

height for both Geldart classes.

In contrast, the normalized particle height is underestimated

by the heterogeneous drag models. The structure-based drag

law by Wang et al. reduced the median for both Geldart types

to a similar height, which yields in a large error for Geldart

B and good results for Geldart D particles. However, the

distribution of particle height is much larger when using the

drag law by Wang et al. compared to the unscaled benchmark

simulation.

123

Computational Particle Mechanics

Table 2 Deviation of the median of the normalized particle heights compared to the benchmark cases at the high gas inlet velocity (for Geldart

class D: uin =5.2umf and for class B: uin =14.2umf) in the wide system Np=100

Drag model Gidaspow [47] (%) Wang et al. [26] (%) EMMS/Matrix [27,28] (%) EMMS/s-s [29,30] (%) Radl et al. [31] (%) Lu et al. [35](%)

Geldart D 12.9−4.7−5.6−6.6−12.9−5.2

Geldart B 5.1−33.0−50.8−38.3−27.3−1.5

The EMMS/Matrix law reduces the normalized particle

height for both particle types. The prediction of the median

for Geldart class D particles is similar to the benchmark

case. The prediction for Geldart class B particles is much

lower and reaches only twice the initial bed height. Look-

ing at the EMMS/steady-state model, for Geldart type B, the

bed height is reduced from four to two times the initial bed

height. For Geldart type D, the effect is less pronounced.

However, less outliers on the top are observed compared to

the EMMS/Matrix model, which is less in agreement with the

benchmark case. Depending on the Geldart class, the model

uses two different correlations. The correlation for Geldart

D type represents the particle dynamics better in this case.

The regression-based Radl model also underestimates the

results for both Geldart classes. The distribution for Gel-

dart D particles gets wider, whereas for Geldart B particles,

more outliers above the upper whisker are observed. The

regression-based model by Lu et al. represents the results

well, but the distribution becomes slightly wider. The wider

distribution of particles in the bed can be observed for the

homogeneous Gidaspow drag law and the heterogeneous

drag laws by Radl et al., Lu et al. and Wang et al., which can

be the statistical influence because the total number of par-

ticles is decreased significantly. For the EMMS models with

Geldart class B particles, the movement of particles is gen-

erally much lower as the models generate the strongest drag

reduction. The results are also summarized in Table 2, where

the relative deviation of the median value of the normalized

particle height compared to the unscaled benchmark simula-

tion as reference is given. The relative deviation is defined as

x=(x−xref)/xref, where xrefers to the analysis variable,

in this case, the median of the normalized bed height.

The effect of the drag models on the fluidization of the par-

ticles is also exemplarily depicted in Fig.8, where snapshots

of the particle distribution at a simulation time of 10 s are

presented for Geldart class D. The simulation with the Wang

law shows only a slight elevation of the particles. The EMMS

group has a similar tendency as well as the filtered drag

laws. Moreover, the time-averaged axial void fraction pro-

files are shown in Fig.9. The void fraction profiles for Geldart

class D particles are given in Fig. 9a. Here, it is also visible

that the heterogeneous drag models all reduce drag a similar

amount, which yields in a slight underprediction of the bed

height, whereas the homogeneous drag model by Gidaspow

leads to overprediction. For Geldart B, shown in Fig.9b, the

Fig. 8 Snapshots taken at 10 s of Geldart class D simulations with

an inlet velocity of 5.2umf in the wide system of Np=100 for the

applied drag models at coarse grain factor of l=10, ashows reference

results for the unscaled benchmark case, bfor the homogeneous drag

law by Gidaspow, cfor the EMMS/steady-state model and dfor the

EMMS/Matrix model. The snapshots for regression-based models by

Radl, Lu, and the structure-based law by Wang are shown in (e), (f),

and (g), respectively

drag reduction for the different heterogeneous drag models

is different and the EMMS/Matrix drag law has the strongest

deviation. The model by Lu et al. predicts the benchmark

case the closest, but still shows slightly different behavior

with lower void fraction below a relative height of 0.4 and

larger void fraction above that height. From the snapshots and

void fraction profiles, the EMMS/Matrix drag law and drag

law by Lu et al. are the most promising for Geldart class D, Lu

et al. is the most promising for Geldart class B, whereas the

others strongly underestimate the normalized particle height

for Geldart class B particles.

4.2.2 Effect of fluidization regime

In fluidized beds, the fluidization regime plays an important

role. While in the bubbling regime, gas bubbles are formed

regularly, accompanied by large fluctuations in pressure and

bed height. Higher slip velocities are reached in the slug-

ging regime, and bubbles transform into large voids, leading

to even larger pressure fluctuations [1,56]. In fluidized beds,

the transition from one fluidization regime to another is often

non-smooth and difficult to determine before conducting

experiments or simulations. Thus, the heterogeneous drag

model must apply to various regimes.

In Fig.10, the normalized particle heights are shown for

the sub-grid drag models and Gidaspow drag law, here exem-

123

Computational Particle Mechanics

(a) Geldart D (b) Geldart B

Fig. 9 Time-averaged void fraction profiles along the relative system

height for Geldart class D and B particles in the wide system (Np=100)

at the large inlet velocity (D: uin =5.2umf,classB:uin =14.2umf)for

a coarse grain factor of l=10. The unscaled benchmark is given as a

reference

Fig. 10 Normalized particle height for Geldart B type particles at the

two fluidization regimes for a coarse grain factor of l=10 in the wide

system Np=100

plary shown for the Geldart B type particle using a coarse

grain factor of 10. The benchmark results are also depicted

for comparison. The results of the simulations for the high

inlet velocity of 14umf show a large underestimation of nor-

malized particle height except for the drag law by Lu et al.,

as already discussed in Sect. 4.2.1. Looking at the low inlet

velocity of 4.2umf, the Gidaspow model overpredicts the bed

height only slightly, indicating that the homogeneous drag

law is also applicable in low-fluidization regimes. However,

the filtered drag models must incorporate this and only apply

slight modifications to the drag model. The drag law by Wang

et al. leads to a large overprediction of bed height in this case.

The EMMS/Matrix drag law and drag law by Radl et al. show

almost no movement of the particles, a narrow distribution,

with only a slight increase in height compared to the ini-

tial bed height. The EMMS/steady-state model gives a slight

decrease compared to the Gidaspow model at low inlet veloc-

ities, indicating that it is suitable for lower gas velocities for

Geldart class B in the bubbling regime. The results are also

summarized in Table 3, where the deviation of the median

compared to the unscaled benchmark simulations is shown.

4.2.3 Influence of system width

At last, the influence of system width is analyzed. As demon-

strated in Sect. 4.1, there is a dependency on the accuracy in

predicting normalized particle height on N.Atlowerratios,a

larger deviation to the unscaled simulation is found. The sub-

grid drag models must reduce the drag accordingly within the

narrower systems.

Figure11 shows the normalized particle heights for the three

system widths with Geldart type B particles at the large inlet

velocity for a coarse grain factor of l=8. Note that l=8is

the maximum coarse grain factor to ensure a minimum of 3

CFD cells in radial direction in the narrow system Np=50.

For the benchmark cases, the behavior in the three systems

shows no significant effect. For l=8 and the Gidaspow

model, the median and the distribution of average particle

height increases from Np=100 to Np=75. From Np=75

to Np=50, the median stays constant. However, the distri-

bution of particle height is much wider for Np=50, reaching

outliers that are larger than nine times the initial bed height.

This indicates different behavior of the particle movement in

the narrow system with large coarse grain factors.

The drag law by Wang et al. reduces the drag strongly

for the system width of Np=100 and Np=75. For

the narrow width of Np=50, the median of the nor-

malized particle height is closer to the benchmark system,

but a wider distribution of the bed height is observed. The

median of the normalized particle height is reduced with

both EMMS models for the narrower system width. How-

ever, bed height is largely underestimated for both laws. The

EMMS/Matrix model predicts a similar normalized particle

height for Np=75 and Np=50, indicating that drag reduc-

tion is stronger in the narrower system. The model by Lu et

al. also gives good agreement for the systems with a sys-

tem diameter of Np=75 and Np=100. All drag models

fail to depict the behavior in the narrow system. Explana-

123

Computational Particle Mechanics

Table 3 Deviation of the median of the normalized particle heights compared to the benchmark cases at the two gas inlet velocities for Geldart

class B in the wide system Np=100

Drag model Gidaspow [47] (%) Wang et al. [26] (%) EMMS/Matrix [27,28] (%) EMMS/s-s [29,30] (%) Radl et al. [31] (%) Lu et al. [35](%)

uin =4.2umf 4.6 18.6−22.53.3−10.46.6

uin =14.2umf 5.1 −33.0−50.8−38.8−27.31 −1.51

Fig. 11 Normalized particle height for Geldart B type particles at the

three different system widths for l=8 at the large inlet velocity (uin =

14umf)

tions on why the narrow system could act differently were

given in Sect. 4.1, and they can also be the reason that the

behavior cannot be corrected with sub-grid drag laws and

other methods might be more applicable, such as interpola-

tion methods to reduce the cell-to-particle ratio. The results

are summarized in Table 4, where the deviation of the median

compared to the unscaled benchmark simulations is given.

4.2.4 Section summary

The effect of the influence of Geldart type, fluidization

regime and system width on the normalized particle height

demonstrates clearly the different behavior of the sub-grid

drag laws applied here, which for most cases for Geldart

class B yields in a strong underprediction of normalized par-

ticle height. For the analyzed cases, the filtered drag model by

Lu et al. performed well except for the narrow system width

Np=50, where the deviation could not be corrected. The

drag law by Wang et al. gave much wider distributions of nor-

malized bed height. However, it could reduce the deviation

for the median looking at the Geldart class D example. Also,

the EMMS/Matrix drag law predicted the Geldart D class

normalized particle height well. For the low inlet velocity,

a drag correction is not necessary. The EMMS/steady-state

model is developed for continuous fluidized bed (CFB) and

bubbling regime, where it performed well at the low inlet

velocity. Similar reasons explain the deviations from the drag

laws by Radl et al. and Wang et al. Moreover, the model by

Wang et al. does not incorporate macroscale markers, e.g.,

gas inlet velocity or system size. The underprediction can

also be explained, when comparing the length scales of the

simulation systems, the models were established or verified

with. Table 7in Sect. A.3 summarizes the ratio of the CFD

cell length to the characteristic cluster length scale Lc.In

the simulation system used here, the ratio reaches maximum

values of 1.36 for Geldart class B and 0.48 for Geldart class

D, while for the other studies, the ratios are larger. This is an

indication that the models are more applicable when larger

length scale ratios are used, such as for larger simulation

systems or for smaller particle types.

4.3 Influence of coarse grain number

In this subsection, the drag law by Lu et al. is investigated on

its accuracy with varying coarse grain factors of l=6,8,10

for Geldart class B and the drag laws by Lu et al., Wang et

al., and the EMMS/Matrix drag law are investigated for Gel-

dart class D and compared to the Gidaspow drag law. As an

example, the deviations of the median and standard devia-

tion of the normalized particle height in the system width of

Np=100 are analyzed. The wider system is chosen to avoid

the error that is caused by the narrower system widths and

not the reduced ratio N.

In Fig.12a and b, the results for the deviation of the median

particle height for both Geldart types are shown. For Geldart

type D, the deviation for the drag law by Gidaspow increases

with increasing coarse grain number. The heterogeneous drag

laws show the opposite trend, leading to a smaller deviation

with a larger coarse grain number, which results in reduced

deviation for coarse grain factors of l=8 and l=10. For

Geldart type B, the trends are similar.

When looking at the deviation of the standard devia-

tion of the particle height signal shown in Fig.12c, the

results for the homogeneous Gidaspow drag law increase

until l=8 and then decrease, similar behavior is observed

for EMMS/Matrix model and the drag law by Lu et al. The

drag law by Wang shows a large deviation for all cases for the

standard deviation of particle height. For the Geldart B par-

ticles shown in Fig.12d, the deviations are increasing with

increasing coarse grain number for the homogeneous model.

The drag law by Lu et al. shows a reduced error except for

l=8.

Looking at the deviations of the median and standard devi-

ation, the results emphasize that the behavior of the drag

models also depends on the coarse grain factor. For a coarse

123

Computational Particle Mechanics

Table 4 Deviation of the median of the normalized particle heights compared to the benchmark cases in the three system width for Geldart class

B at the large inlet velocity of uin =14.2umf for l=8

Drag model Gidaspow [47] (%) Wang et al. [26] (%) EMMS/Matrix [27,28] (%) EMMS/s-s [29,30] (%) Radl et al. [31] (%) Lu et al. [35](%)

Np=50 11.6−11.7 −47.7 −27.4 −3.315.2

Np=75 10.0−34.5 −49.9 −33.4 −19.9−2.3

Np=100 4.1−37.6 −52.9 −38.1 −26.8−1.1

Fig. 12 Absolute values of the

relative deviation of the median

and standard deviation of

particle height to the unscaled

Benchmark case for Geldart

class B particles (a,b)and

Geldart class D (c,d)atthelarge

gas inlet velocity (for Geldart

class D: uin =5.2umf and for

class B: uin =14.2umf)for

different coarse grain factors lin

the wide system of Np=100

(a) Geldart D: Median (b) Geldart B: Median

grain factor of l=6, the error for using the Gidaspow drag

law is generally low, so applying sub-grid drag models, in

this case, leads to a larger error. However, for l=8 and

l=10, the error can be reduced with some sub-grid laws.

5 Conclusions

In this work, different heterogeneous drag models were

analyzed toward their applicability to incorporate sub-grid

behavior. Firstly, the results for the drag law by Gidaspow

were analyzed. Different behavior of Geldart class D and B

was observed, when predicting the normalized particle height

with decreasing ratio N. The deviation to the benchmark case

is larger at the large gas inlet velocity for both Geldart types,

but more pronounced with the class D particles. The increase

of particle height depends not only on the ratio N, but also on

the absolute value of system width as the investigated narrow

width (Np=50) showed a larger deviation at similar ratios

For the simulations with large coarse grain factors of l=

6−10 and subsequently low ratios N, sub-grid drag laws were

applied. Their performances in predicting particle height

using configurations of different particles type, gas inlet

velocity, and system width were analyzed. Looking at the two

Geldart types, the simulations with the sub-grid drag models

obtained different results, in which for the Geldart class D the

error was less pronounced and for class B particles, the nor-

malized bed height was strongly underestimated except for

the sub-grid law by Lu et al. Looking at the two inlet veloc-

123

Computational Particle Mechanics

ities, a lower error at the low inlet velocity of uin =4.2umf,

except for the drag law by Lu et al. was observed compared to

results at the large inlet velocity of uin =14.2umf.Thesim-

ulation results obtained with the narrowest system were not

improved by any drag law, either due to a limiting ratio Nor

due to an insufficient number of CFD cells. The latter can be

improved by using interpolation mapping methods to calcu-

late void fractions and slip velocities, which will need further

evaluation. Lastly, the results of the median and standard

deviation of the normalized particle height were analyzed

for the sub-grid drag laws with the lowest deviations for dif-

ferent coarse grain factors, which emphasized the different

behavior of the sub-grid drag laws with different coarse grain

factors.

The predictions not only vary depending on the Geldart type,

inlet velocity, and system width, but also on the coarse grain

factor. If a drag law is chosen, that does not fit the system’s

configurations, it can lead to an overcorrection of the problem

and might increase the error. This work also shows that there

is not one drag law that can be applied for all fluid regimes,

Geldart types, and system configurations. Depending on the

type of simulation, the sub-grid drag laws are a useful tool

within the CG-DEM–CFD approach but always need some

evaluation and comparison to other drag models. As men-

tioned by [57], the physical phenomena causing the micro-,

meso-, and macroscale effects in fluidized beds are not fully

understood yet, so a general drag formulation needs far more

research.

Acknowledgements Financial support was funded by the Deutsche

Forschungsgemeinschaft (DFG, German Research Foundation): Project

ID 456827728. We would like to express our gratitude to Prof. Wei

Wang from the Institute of Process Engineering, Chinese Academy of

Sciences for providing the EMMS Software and correlations for the

EMMS/steady-state model used in this study.

Funding Open Access funding enabled and organized by Projekt

DEAL.

Declarations

Conflict of interest On behalf of all authors, the corresponding author

states that there is no conflict of interest.

Open Access This article is licensed under a Creative Commons

Attribution 4.0 International License, which permits use, sharing, adap-

tation, distribution and reproduction in any medium or format, as

long as you give appropriate credit to the original author(s) and the

source, provide a link to the Creative Commons licence, and indi-

cate if changes were made. The images or other third party material

in this article are included in the article’s Creative Commons licence,

unless indicated otherwise in a credit line to the material. If material

is not included in the article’s Creative Commons licence and your

intended use is not permitted by statutory regulation or exceeds the

permitted use, you will need to obtain permission directly from the copy-

right holder. To view a copy of this licence, visit http://creativecomm

ons.org/licenses/by/4.0/.

A Appendix

A.1 Governing equations

All simulations are solved with the open-source software

MFiX and the model description is based on the implemented

equations, which can be found in [42]. The solid phase is

described with DEM, position, direction, and velocity of each

particle are calculated with the following equations:

dvp(t)

dt=mpg+Fh

p(t)+Fc

p(t)(13)

dωp(t)

dt=T,(14)

where Fh

p(t)is the total hydrodynamic force acting on

the particle, and Fc

p(t)is the contact force resulting from

particle–particle or wall–particle contacts. Tstands for the

sum of all torques acting on the particle.

In this study, the contact force is calculated using the linear-

spring-dashpot model. The contact force is the sum of a

normal force and a tangential force, and it is calculated for

each particle pair or particle–wall pair ij such as

n,ij =knδijnij −γnvn,ij (15)

t,ij =ktδijtij −γtvt,ij if Fc

t,ij ≤μ|Fc

n,ij|

−μ|Fc

n,ij|tij otherwise (16)

with knas the normal stiffness, ktas the tangential stiff-

ness. δij is the overlap between the particle–particle or

particle–wall pair. nij and tij are the normal and tangential

unit vectors, γnand γtare normal and tangential damping

coefficients. vn,ij and vt,ij is the relative normal and relative

tangential velocities. The hydrodynamic force consists of the

pressure gradient force and the drag force:

p=F∇p

p+Fd

p=−Vp∇pg+βVp

(1−εf)(uf−vp). (17)

The fluid flow is solved with the gas-phase continuity and

momentum conservation as follows:

∂εfρf

∂t+∇·(ε fρfuf)=0 (18)

and

∂

∂t(ε fρfuf)+ρf∇·(ε fufuf)

=εf∇pf+∇·(ε fτf)+εfρfg−Ifp (19)

with pfas the gas-phase pressure, τfas the gas-phase stress

tensor, and Ifp as the momentum transfer term between the

gas and the particles or solid phase p.

123

Computational Particle Mechanics

A.2 Equations for sub-grid drag models

Drag law by Wang et al. [26]

The drag occurring in the dense phase is calculated by using

the homogeneous drag model with the dense phase parame-

ters:

βd=β(εf,d,εp,d,Red)(20)

with εf,das the void fraction and εp,dthe solid volume frac-

tion in the dense phase. The Reynolds number is defined as:

Red=ρfdp|Ud|

μf

(21)

with Udas the superficial gas velocity of the dense phase,

the is void fraction defined as:

εf,d=1−εp,d.(22)

With these, the drag force within the emulsion phase can be

calculated with:

Fd=βd

εf,d

.(23)

For any computational grid cell, a varying part of it is occu-

pied by the dense phase, and the drag force within the cell

amounts to:

Fcell =βd

εg,d

εp

εp,d

,(24)

where εp

εp,dis the emulsion phase fraction in the cell.

Now, only εp,d(εf,d=1−εp,d) and Udneed to be deter-

mined. They are derived from experimental correlations. εp,d

can be calculated with the following equations:

εp,d

1−εmf =1−0.14Re0.4Ar−0.13 (25)

with

Re =ρfdp||uf−vp|εg−umf|

μf

(26)

and

Ar =ρfd3

p(ρp−ρf)g

μ2

,(27)

where εmf is the voidage at minimum fluidization. With

εf,d,Udcan be calculated with the following relation:

|Ud|

vt=εn

g,d(28)

with vtbeing the particle terminal velocity and with the expo-

nent nbeing:

n=ln(umf/vt)

ln(εmf).(29)

The drag coefficient for each computational cell can now

be calculated with:

βcell =|Fcell|

(uf−vp)(30)

EMMS/Matrix [27,28]

The EMMS/Matrix heterogeneity index is read into the soft-

ware by bilinear interpolation. The heterogeneity index as a

function of void fraction and slip velocity is given in Fig.13.

EMMS/steady-state [29,30]

In the following, the equations used for the EMMS/steady-

state model are stated. Equation31 is used for simulations

with Geldart D particles at an inlet velocity of uin =2.5umf,

Eq.32 for simulations with Geldart D particles at an inlet

velocity of uin =5.2umf. Equation 33 is used for simulations

with Geldart B particles at an inlet velocity of uin =4.2umf,

Eq.34 for simulations with Geldart B particles at an inlet

velocity of uin =14.2umf.

β=⎧

⎪

⎨

⎪

⎩

1,(1−εf)<0.4237

exp(y), 0.4237 <(1−εf)<05496

1,(1−εf)>0.5496

(31)

with

y=−153498(1−εf)5+387781.3(1−εf)4

−388814.1−εf)3

+193606.5(1−εf)2−47922.4(1−εf)+4721.26).

β=⎧

⎪

⎨

⎪

⎩

1,(1−εf)<0.2707

exp(y), 0.2707 <=(1−εf)<0.5467

1,(1−εf)>0.5467

(32)

with

y=−3552961(1−εf)7+9904061(1−εf)6

+−11714213.89(1−εf)5

123

Computational Particle Mechanics

Fig. 13 EMMS/Matrix contour

plot for heterogeneity index as a

function of slip velocity and

void fraction for Geldart B at the

low inlet velocity uin =4.2umf

(a) and at the high inlet velocity

uin =14.2umf (b)aswellasfor

Geldart D at the low inlet

velocity uin =2.5umf (c)and

the high inlet velocity

uin =5.2umf (d)

(a) Geldart B (b) Geldart B

+7621751(1−εf)4+−2946519(1−εf)3

+676886.6(1−εf)2

+−85570.2(1−εf)+4593.392

β=⎧

⎪

⎨

⎪

⎩

1,(1−εf)<0.4027

exp(y), 0.4027 <=(1−εf)<0.5616

1,(1−εf)>0.5616

(33)

with

y=−4606259(1−εf)6+13304060(1−εf)5

−15956281.75(1−εf)4

+10172982(1−εf)3−3636637(1−εf)2

+691187.6(1−εf)−54569.5.

β=⎧

⎪

⎨

⎪

⎩

1,(1−εf)<0.1388

exp(y), 0.1388 <=(1−εf)<0.5796

1,(1−εf)>0.5796

(34)

with

y=−319008(1−εf)7+762909.6(1−εf)6

−755893(1−εf)5

+401771.9,(1−εf)4−123542(1−εf)3

+21938.22(1−εf)2

−2085.56(1−εf)+81.89.

Drag law by Radl et al. [31]

βis calculated as:

β=1−ff

, εph(εp)ccorr(l, εp)(35)

with

ff

, εp=1

a(εp)(Lc/ f)+1(36)

and

ccorr =exp(0.05(l−1)) (37)

with a(εp)and h(εp)being piecewise continuous alge-

braic expressions presented in Tables 5and 6. The function

ccorr was introduced to account for particle coarsening and

has been validated for coarse grain factors up to l=8.

The characteristic cluster length scale Lcis given by:

Lc=v2

gFrn

p(38)

with vtbeing the terminal particle velocity and Fr being the

particle Froude number. The exponent nwas determined to

be −2

3in order to minimize the influence of the filter size on

the simulated parameters.

123

Computational Particle Mechanics

Table 5 Piecewise continuous

algebraic expression for

a(εp)=

a0,m+a1,mεp−εp,a,m+

a2,mεp−εp,a,m2+

a3,mεp−εp,a,m3

mRange a0,ma1,ma2,ma3,mεp,a,m

1εp<0.016 21.51 0 0 0 0

20.016 <εp<0.10 1.96 29.40 164.91 −1923 0

30.10 <εp<0.18 4.63 4.68 −412.04 2254 0.10

40.18 <εp<0.25 3.52 −17.99 128.80 −603 0.18

50.25 <εp<0.40 2.68 −8.82 2.18 112.33 0.25

60.40 <εp1.79 0 0 0 0

Table 6 Piecewise continuous

algebraic expression for

h(εp)=

h0,m+h1,mεp−εp,a,m+

h2,mεp−εp,a,m2+

h3,mεp−εp,a,m3

mRange h0,mh1,mh2,mh3,mεp,a,m

1εp<0.03 0 7.97000

20.03 <εp<0.08 0.239 4.64 −4.41 253.63 0.03

30.08 <εp<0.12 0.492 6.10 33.63 −789.60.08

40.12 <εp<0.18 0.739 5.01 −61.1 310.80.12

50.18 <εp<0.34 0.887 1.03 −5.17 5.99 0.18

60.34 <εp<0.48 0.943 −0.170 −2.29 −9.12 0.34

70.48 <εp<0.55 0.850 −1.35 −6.13 −132.60.48

80.55 <εp<0.60 0.680 −2.34 −225.20 0.55

Table 7 CFD cell size and characteristic cluster length scale of the par-

ticle systems in literature. These simulations were used to verify the

applied drag models. For comparison, the length scales for the system

in this study are also shown for Geldart B and Geldart D particles for a

coarse grain factor of l=10. The terminal velocity is calculated with

Eq.46, the characteristic cluster length scale with Eq. 38

Model Simulation model Verification vt/(m/s) dx/m Lc/m dx/Lc

Wang et al. [26] TFM Experimental data based on Werther & Wein 4.53 5.00 ×10−24.94 ×10−310.12

Wang et al. [26] TFM Experimental data based on Johnson et al. 24.68 2.26 ×10−12.69 ×10−28.42

Radl et al. [31] CGDEM–CFD (l=6) DEM–CFD 0.26 2.00 ×10−33.34 ×10−45.98

Lu et al. [35] CGDEM–CFD (dx =6dp) DEM–CFD 0.26 4.50 ×10−43.34 ×10−41.35

Lu et al. [35] CGDEM–CFD (dx =12dp) DEM–CFD 0.26 9.00 ×10−43.34 ×10−42.69

Lu et al. [35] CGDEM–CFD (dx =15dp) DEM–CFD 0.26 1.13 ×10−33.34 ×10−43.37

EMMS/Matrix [27]TFM 0.08 9.69 ×10−31.26 ×10−476.96

Simulation case

Geldart B CGDEM–CFD (l=10) DEM–CFD 11.35 2 ×10−21.47 ×10−21.36

Geldart D CGDEM–CFD (l=10) DEM–CFD 75.66 4 ×10−28.29 ×10−20.48

Drag law by Lu et al. [35]

Hdis calculated with:

Hd=2−1.99(1−e−α(u∗

slip−u0)P+X,u∗

slip >u0

1,u∗

slip ≤u0,εp>0.59,εp<0.01

(39)

with the parameters α, u0,P,X,fand u∗

slip as the filtered,

dimensionless slip velocity:

α=(a1+a2εp+a3ε2

p+a4ε3

p+a5ε4

p)(1−ea19εp)

1+ea20(εp−0.55)

1+a6

f∗+a7

( f∗)21+a8

(u∗

slip)2(40)

u0=a9+a10εp

0.01 +εa11

p1+a12

f∗+a13

( f∗)2(41)

P=a14 +a15εp+a16ε2

p1+a17

f∗+a18

( f∗)2(42)

X=a21 +a22 +a23f∗+a24εp+a25u∗

slipUf

(43)

f∗=2gf

(44)

u∗

slip =uslip

(45)

123

Computational Particle Mechanics

with the filtered slip velocity uslip =|uf−vp|and the ter-

minal particle velocity vt:

vt=gd2

p(ρg−ρp)

18 ·μg

(46)

and fas the filter size:

f=2·V1/3

cell (47)

where Vcell is the volume of the fluid cell. The regressed

coefficients are:

a1−25 =[4.802706,12.753237,67.982539,

96.535179,104.148042,0.014874,0.200902,

0.00102,0.038998,1.492272,

0.317903,0.121642,0.024993,0.155943,

4.673687,18.12709,0.016406,0.248029,

707.455212,5.49012,0.086153,0.614471,

0.174122,0.546865,0.170773].

A.3 Comparison of length scales in the models

The applied models in this study were all verified based on

simulations of fluidized beds and compared either to experi-

mental or fine grid simulation data. For the coarse simulations

in these studies, the CFD cell length dx and the character-

istic cluster length scale Lcfrom Eq. 38 are summarized in

Table 7. For comparison, the length scales used in this study

for both particle classes (Geldart class B and class D) for a

coarse grain factor of l=10 are also shown.

References

1. Grace J, Bi X, Ellis N (2020) Essentials of fluidization technology.

Wiley, New York. https://doi.org/10.1002/9783527699483

2. Di Renzo A, Scala F, Heinrich S (2021) Recent advances in

fluidized bed hydrodynamics and transport phenomena-progress

and understanding. Processes 9(4):639. https://doi.org/10.3390/

pr9040639

3. Chew JW, LaMarche WCQ, Cocco RA (2022) 100 years of scaling

up fluidized bed and circulating fluidized bed reactors. Pow-

der Technol 409:117813. https://doi.org/10.1016/j.powtec.2022.

117813

4. Vollmari K, Oschmann T, Kruggel-Emden H (2017) Mixing quality

in mono- and bidisperse systems under the influence of parti-

cle shape: a numerical and experimental study. Powder Technol

308:101–113. https://doi.org/10.1016/j.powtec.2016.11.072

5. Alobaid F, Almohammed N, Massoudi Farid M, May J, Rößger

P, Richter A, Epple B (2022) Progress in CFD simulations of flu-

idized beds for chemical and energy process engineering. Prog

Energy Combust Sci 91:100930. https://doi.org/10.1016/j.pecs.

2021.100930

6. Kruggel-Emden H, Kravets B, Suryanarayana MK, Jasevicius R

(2016) Direct numerical simulation of coupled fluid flow and

heat transfer for single particles and particle packings by a LBM-

approach. Powder Technol 294:236–251. https://doi.org/10.1016/

j.powtec.2016.02.038

7. Deen NG, Kriebitzsch SHL, Hoef MA, Kuipers JAM (2012) Direct

numerical simulation of flow and heat transfer in dense fluid-

particle systems. Chem Eng Sci 81:329–344. https://doi.org/10.

1016/j.ces.2012.06.055

8. Lu L, Liu X, Li T, Wang L, Ge W, Benyahia S (2017) Assessing the

capability of continuum and discrete particle methods to simulate

gas-solids flow using DNS predictions as a benchmark. Powder

Technol 321:301–309. https://doi.org/10.1016/j.powtec.2017.08.

034

9. Esteghamatian A, Bernard M, Lance M, Hammouti A, Wachs A

(2017-06) Micro/meso simulation of a fluidized bed in a homoge-

neous bubbling regime. Int J Multiph Flow 92:93–111. https://doi.

org/10.1016/j.ijmultiphaseflow.2017.03.002

10. Ding J, Gidaspow D (1990) A bubbling fluidization model using

kinetic theory of granular flow. AIChE J 36(4):523–538. https://

doi.org/10.1002/aic.690360404

11. Patil DJ, Sint Annaland M, Kuipers JAM (2005) Critical compar-

ison of hydrodynamic models for gas-solid fluidized beds-part II:

freely bubbling gas-solid fluidized beds. Chem Eng Sci 60(1):73–

84. https://doi.org/10.1016/j.ces.2004.07.058

12. Patil DJ, Sint Annaland M, Kuipers JAM (2005) Critical compar-

ison of hydrodynamic models for gas-solid fluidized beds-part I:

bubbling gas–solid fluidized beds operated with a jet. Chem Eng

Sci 60(1):57–72. https://doi.org/10.1016/j.ces.2004.07.059

13. Yang L, Padding JTJ, Kuipers JAMH (2017) Investigation of col-

lisional parameters for rough spheres in fluidized beds. Powder

Technol 316:256–264. https://doi.org/10.1016/j.powtec.2016.12.

090

14. Hoomans BPB (2000) Granular dynamics of gas–solid two-phase

flows. Ph.D. thesis, Universiteit Twente, The Netherlands

15. Schneiderbauer S, Puttinger S, Pirker S (2013-11) Comparative

analysis of subgrid drag modifications for dense gas-particle flows

in bubbling fluidized beds. AIChE J 59(11):4077–4099. https://doi.

org/10.1002/aic.14155

16. Di Renzo A, Napolitano ES, Di Maio FP (2021) Coarse-grain dem

modelling in fluidized bed simulation: a review. Processes 9(2):1–

30. https://doi.org/10.3390/pr9020279

17. Lin J, Luo K, Wang S, Hu C, Fan J (2020) An augmented coarse-

grained CFD-DEM approach for simulation of fluidized beds. Adv

Powder Technol 31(10):4420–4427. https://doi.org/10.1016/j.apt.

2020.09.014

18. Stroh A, Daikeler A, Nikku M, May J, Alobaid F, Bohnstein MV,

Ströhe J, Epple B, Stroh A, Daikeler A, Nikku M, May J, Alobaid

F (2018) Coarse grain 3D CFD-DEM simulation and validation

with capacitance probe measurements in a circulating fluidized

bed. Chem Eng Sci. https://doi.org/10.1016/j.ces.2018.11.052

19. Napolitano ES, Di Renzo A, Di Maio FP (2022) Coarse-grain

DEM-CFD modelling of dense particle flow in gas–solid cyclone.

Sep Purif Technol 287:120591. https://doi.org/10.1016/j.seppur.

2022.120591

20. Scott L, Borissova A, Di Renzo A, Ghadiri M (2022) Applica-

tion of coarse-graining for large scale simulation of fluid and

particle motion in spiral jet mill by CFD-DEM. Powder Technol

411:117962. https://doi.org/10.1016/j.powtec.2022.117962

21. Igci Y, Iv ATA, Sundaresan S, Brien TO (2008) Filtered two-fluid

models for fluidized gas-particle suspensions. AIChE J. https://doi.

org/10.1002/aic.11481

22. Agrawal K, Loezos PN, Symlal M, Sundaresan S (2001) The role

of meso-scale structures in rapid gas–solid flows. J Fluid Mech

445:151–185. https://doi.org/10.1017/S0022112001005663

23. Zhang DZ, VanderHeyden WB (2002) The effects of mesoscale

structures on the macroscopic momentum equations for two-phase

123

Computational Particle Mechanics

flows. Int J Multiph Flow 28(5):805–822. https://doi.org/10.1016/

S0301-9322(02)00005-8

24. Wang J (2009) A review of Eulerian simulation of Geldart a parti-

cles in gas-fluidized beds. Ind Eng Chem Res 48(12):5567–5577.

https://doi.org/10.1021/ie900247t

25. O’Brien TJ, Syamlal M (1993) Particle cluster effects in the numeri-

cal simulation of a circulating fluidized bed. Circ Fluid Bed Technol

IV:367–372

26. Wang J, van der Hoef MA, Kuipers JAM (2010) Coarse grid simu-

lation of bed expansion characteristics of industrial-scale gas - solid

bubbling fluidized beds. Chem Eng Sci 65(6):2125–2131. https://

doi.org/10.1016/j.ces.2009.12.004

27. Wang W, Li J (2007) Simulation of gas-solid two-phase flow by

a multi-scale CFD approach-of the EMMS model to the sub-grid

level. Chem Eng Sci 62(1):208–231. https://doi.org/10.1016/j.ces.

2006.08.017

28. Ge W, Li J (2002) Physical mapping of fluidization regimes-the

EMMS approach. Chem Eng Sci 57(18):3993–4004. https://doi.

org/10.1016/S0009-2509(02)00234-8

29. Tian Y, Lu B, Li F, Wang W (2020) A steady-state EMMS drag

model for fluidized beds. Chem Eng Sci 219:115616. https://doi.

org/10.1016/j.ces.2020.115616

30. Lu B, Zhang N, Wang W, Li J (2012) Extending EMMS-based

models to CFB boiler applications. Particuology 10(6):663–671.

https://doi.org/10.1016/j.partic.2012.06.003

31. Radl S, Sundaresan S (2014) A drag model for filtered Euler–

Lagrange simulations of clustered gas-particle suspensions. Chem

Eng Sci 117:416–425. https://doi.org/10.1016/j.ces.2014.07.011

32. Igci Y, Iv ATA, Sundaresan S, Brien TO (2008) Filtered two-fluid

models for fluidized gas-particle suspensions. AIChE J. https://doi.

org/10.1002/aic

33. Lu B, Wang W, Li J (2009) Searching for a mesh-independent sub-

grid model for CFD simulation of gas-solid riser flows. Chem Eng

Sci 64(15):3437–3447. https://doi.org/10.1016/j.ces.2009.04.024

34. Jiang Y, Kolehmainen J, Gu Y, Kevrekidis YG, Ozel A, Sun-

daresan S (2019) Neural-network-based filtered drag model for

gas-particle flows. Powder Technol 346:403–413. https://doi.org/

10.1016/j.powtec.2018.11.092

35. Lu L, Gao X, Dietiker J-F, Shahnam M, Rogers WA (2021) Devel-

opment of a filtered CFD-DEM drag model with multiscale markers

using an artificial neural network and nonlinear regression. Ind Eng

Chem Res. https://doi.org/10.1021/acs.iecr.1c03644

36. Du C, Han C, Yang Z, Wu H, Luo H, Niedzwiecki L, Lu B, Wang

W (2022) Multiscale CFD simulation of an industrial diameter-

transformed fluidized bed reactor with artificial neural network

analysis of EMMS drag markers. Ind Eng Chem Res 61(24):8566–

8580. https://doi.org/10.1021/acs.iecr.2c00396

37. Yang Z, Lu B, Wang W (2021) Coupling artificial neural network

with EMMS drag for simulation of dense fluidized beds. Chem Eng

Sci 246:117003. https://doi.org/10.1016/j.ces.2021.117003

38. Lu L, Xu J, Ge W, Yue Y, Liu X, Li J (2014) EMMS-based discrete

particle method (EMMS-DPM) for simulation of gas–solid flows.

Chem Eng Sci 120:67–87. https://doi.org/10.1016/j.ces.2014.08.

004

39. Jurtz N, Kruggel-Emden H, Cocco R (2020) Impact of contact scal-

ing and drag calculation on the accuracy of coarse-grained discrete

element method. Chem Eng Technol. https://doi.org/10.1002/ceat.

202000055

40. Adnan M, Sun J, Ahmad N, Wei JJ (2020) Multiscale modeling of

bubbling fluidized bed reactors using a hybrid Eulerian–Lagrangian

dense discrete phase approach. Powder Technol 376:296–319.

https://doi.org/10.1016/j.powtec.2020.07.111

41. Wang X, Chen K, Kang T, Ouyang J (2020) A dynamic coarse

grain discrete element method for gas–solid fluidized beds by

considering particle-group crushing and polymerization. Appl

Sci(Switzerland) 10:6. https://doi.org/10.3390/app10061943

42. Syamlal M (1994) MFIX documentation: users manual. Techni-

cal report, EG and G Technical Services of West Virginia, Inc.,

Morgantown

43. Syamlal M, et al (1998) Mfix documentation: numerical tech-

nique. National Energy Technology Laboratory, Department

of Energy, Technical Note Nos. DOE/MC31346-5824 and

NTIS/DE98002029 see also http://www.mfix.org

44. Benyahia S, Galvin JE (2010) Estimation of numerical errors

related to some basic assumptions in discrete particle methods.

Ind Eng Chem Res 49(21):10588–10605. https://doi.org/10.1021/

ie100662z

45. Cai R, Zhao Y (2020) An experimentally validated coarse-grain

DEM study of monodisperse granular mixing. Powder Technol

361:99–111. https://doi.org/10.1016/j.powtec.2019.10.023

46. Brandt V, Grabowski J, Jurtz N, Kraume M, Kruggel-Emden H

(2023) A benchmarking study of different DEM coarse graining

strategies. Powder Technol 426:118629. https://doi.org/10.1016/j.

powtec.2023.118629

47. Gidaspow D (1994) Flow and fluidization: continuum and kinetic

theory descriptions. Academic Press, Boston. https://doi.org/10.

1016/C2009-0-21244-X

48. Ergun S (1952) Fluid through packed columns. Chem Eng Prog

48(2):89–94

49. Wen CY, Yu YH (1966) Mechanics of fluidization. Chem Eng Prog

Symp Ser 62(1):100–111

50. Sarkar A, Milioli FE, Ozarkar S, Li T, Sun X, Sundaresan S

(2016) Filtered sub-grid constitutive models for fluidized gas-

particle flows constructed from 3-D simulations. Chem Eng Sci

152:443–456. https://doi.org/10.1016/j.ces.2016.06.023

51. Gao X, Li T, Sarkar A, Lu L, Rogers WA (2018) Development and

validation of an enhanced filtered drag model for simulating gas-

solid fluidization of Geldart a particles in all flow regimes. Chem

Eng Sci 184:33–51. https://doi.org/10.1016/j.ces.2018.03.038

52. Kruggel-Emden H, Stepanek F, Munjiza A (2010) A study on

adjusted contact force laws for accelerated large scale discrete ele-

ment simulations. Particuology 8(2):161–175. https://doi.org/10.

1016/j.partic.2009.07.006

53. Clarke DA, Sederman AJ, Gladden LF, Holland DJ (2018) Investi-

gation of void fraction schemes for use with CFD-DEM simulations

of fluidized beds. Ind Eng Chem Res 57(8):3002–3013. https://doi.

org/10.1021/acs.iecr.7b04638

54. Garg R, Galvin J, Li T, Pannala S(2012) Documentation of

open-source mfix–dem software for gas–solids flows. https://mfix.

netl.doe.gov/download/mfix/mfix_current_documentation/dem_

doc_2012-1.pdf

55. Goldschmidt MJV, Link JM, Mellema S, Kuipers JAM (2003) Dig-

ital image analysis measurements of bed expansion and segregation

dynamics in dense gas-fluidised beds. Powder Technol 138(2):135–

159. https://doi.org/10.1016/j.powtec.2003.09.003

56. Thornton C, Yang F, Seville J (2015) A DEM investigation

of transitional behaviour in gas-fluidised beds. Powder Technol

270:128–134. https://doi.org/10.1016/j.powtec.2014.10.017

57. Wang W, Lu B, Geng J, Li F (2020) Mesoscale drag modeling: a

critical review. Curr Opin Chem Eng 29:96–103. https://doi.org/

10.1016/j.coche.2020.07.001

Publisher’s Note Springer Nature remains neutral with regard to juris-

dictional claims in published maps and institutional affiliations.

123