Coarse-Grain DEM Modelling in Fluidized Bed Simulation: A Review

: In the last decade, a few of the early attempts to bring CFD-DEM of fluidized beds beyond the limits of small, lab-scale units to larger scale systems have become popular. The simulation ca-pabilities of the Discrete Element Method in multiphase flow and fluidized beds have largely bene-fitted by the improvements offered by coarse graining approaches. In fact, the number of real particles that can be simulated increases to the point that pilot-scale and some industrially relevant systems become approachable. Methodologically, coarse graining procedures have been introduced by various groups, resting on different physical backgrounds. The present review collects the most relevant contributions, critically proposing them within a unique, consistent framework for the der-ivations and nomenclature. Scaling for the contact forces, with the linear and Hertz-based approaches, for the hydrodynamic and cohesive forces is illustrated and discussed. The orders of magnitude computational savings are quantified as a function of the coarse graining degree. An over-view of the recent applications in bubbling, spouted beds and circulating fluidized bed reactors is presented. Finally, new scaling, recent extensions and promising future directions are discussed in perspective. In addition to providing a compact compendium of the essential aspects, the review aims at stimulating further efforts in this promising field. for corrections to the conventional drag models for in risers and coarse graining


Introduction
The fluidized bed technology is at the very heart of a very broad number of industrial processes, ranging from chemical transformations in reactors (energy conversion and storage, (bio-)oil refining, chemical synthesis, polymerization) to physical operations (solids mixing or separation, drying, coating, agglomeration) [1].Specialized applications can be found in the development of bioartificial organs [2,3].
Modelling of the gas-and liquid-solid flow in fluidized beds has traditionally relied on the similarity of the fluidized solid motion with a pseudo-liquid phase, through the concept of a gas-solid suspension.Even from a lexical point of view, the terms "fluidization", "bubbling", "emulsion phase", "floating" and "sinking" of objects in the suspended particle bed recall the behavior of liquids.Based on this idea, the Two-Fluid Model (TFM) concept [4] was introduced, complemented by the kinetic theory of granular flow (KTGF), giving rise to groundbreaking success in the simulation of large-scale systems, as summarized in the book by Gidaspow [5].At the other extreme, i.e., the small-scale interaction between the fluid and individual particles, discrete-continuum (or Eulerian-Lagrangian) simulations based on the discrete element method appeared, thanks to the progress in the computational power, giving rise to the CFD-DEM (Computational Fluid Dynamics-Discrete Element Method) approach [6].The recognition of the discrete nature of the granular medium proved to be a key ingredient and allowed important new features of fluidized systems to be captured, as discussed in several reviews (e.g., [7][8][9][10][11][12]).
The most significant drawback of CFD-DEM is the feasibility of realistic simulations, as the combination of practical sizes, number of particles and required time to simulate very quickly saturates the computing capability of the most powerful supercomputer.As examples of the achievable scale in terms of numbers of particles, CFD-DEM simulations reached 4.5 million particles in bubbling beds [13] and up to 25 million particles using the power of the GPU [14].In a recent review on methods for granular and particle-fluid flow, with focus on strategies for engineering applications, Ge et al. [15] noted that the typical time step is 10 −5 -10 −6 s and an engineering system may contain 10 8 -10 12 particles, with characteristic flow or transport time ranges of 10 1 -10 4 s (see also [14]).The necessary time and hardware resources to perform direct DEM simulations and post-process the results of such systems are prohibitively costly and hence they remain infeasible.
In the last years, the push to save CPU time while maintaining the discrete nature of the simulated elements has been a very strong driver for novel strategies and methods.From the physical point of view, Lu et al. [16] noted that in most fluidized bed applications, only the collective behavior of the particles is of primary interest, rather than their individual trajectories.This forms the basis of practical demand for the application of coarse-graining to DEM.The idea behind coarse-grained DEM (or CG-DEM) is to substitute the actual particles by a smaller number of representative elements, whose behavior shall be equivalent in full to the original ones.These representative particles, often called "parcels", become the elements whose motion is tracked by the simulations.They are not required to exist in reality, but their study provides useful data and insight on the real system they represent.The concept of parcels is very common in the context of the MP-PIC method, another Eulerian-Lagrangian technique which falls between TFM and CFD-DEM; the main difference is that solid stresses are treated by indirect models rather than explicitly like in DEM.
Previous attempts to modify the contact parameters to speed up DEM simulations is extremely common for particulate systems in the collisional or dilute regime.At the risk of introducing artifacts, some authors investigated different model laws to speed up the simulation; see e.g., [17].In the coarse graining method, the objective is to scale size and parameters, and adapt models, so that the properties at the scale of the represented particles are kept constant.Among the early ideas of coarse-grained DEM, the following approaches can be mentioned: the imaginary sphere model [18], Similar Particle Assembly [19,20], Representative Particle Model [21]; the work of Patankar and Joseph [22], in which one model of the Lagrange treatment of the solid phase is based on a parcel-like DEM model, deserves to be mentioned.The two CG-DEM models that found widespread adoptions were proposed by Sakai et al. [23] for pneumatic conveying and by Bierwisch et al. [24] in the context of simulation of the cavity filling process.In the same period, Radl et al. [25] proposed scaling rules for parcels in dense gas-particle flows and Hilton and Cleary [26,27] showed an application to fluidized beds.The idea of using scaling rules to save on computational time was also tested [28], with limited success, as with exact scaling the simulation speed-up came at the cost of longer times to simulate, frustrating the ambitions of the method.Other interesting results related to similarity and simulations are found in [29][30][31].Overall, application of the coarse graining method to DEM simulations is attracting a quickly growing attention in recent times, as shown by the yearly number of papers on the subject (see Figure 1).Early applications have been reported in commercial software packages [32] and are available in recent versions, see e.g., [33][34][35].The feature is available also in the last version (20.4) of the open-source software MFiX-DEM [36].
The present contribution aims at critically reviewing the work done so far, presenting the approaches that found widespread adoption within a consistent and coherent framework, occasionally comparing the formulations; quantifying the theoretical computational savings as well as discussing the most critical issues affecting accuracy and reliability; presenting a picture of the state-of-the-art and finally discussing the most recent extensions and the questions still open.
A final introductory note to the reader: One should be aware that the topics covered here do not include other methods incidentally also known as "coarse graining" in DEM, which are used to migrate quantities from the micro-to the macroscale (like e.g., [37]) or smoothing discrete quantities for coupling ( [38]).
Processes 2021, 9, x FOR PEER REVIEW 3 of 30 tional savings as well as discussing the most critical issues affecting accuracy and reliability; presenting a picture of the state-of-the-art and finally discussing the most recent extensions and the questions still open.
A final introductory note to the reader: One should be aware that the topics covered here do not include other methods incidentally also known as "coarse graining" in DEM, which are used to migrate quantities from the micro-to the macroscale (like e.g., [37]) or smoothing discrete quantities for coupling ( [38]).

CFD-DEM Modelling for Fluidized Beds
In CFD-DEM, the equations governing the two phases are solved separately, i.e., sequentially.Since the characteristic times can be different, distinct time-steps are used for the CFD and the DEM parts, so that generally the fluid flow remains constant for a number of DEM (smaller) time-steps.
The fluid phase flow is solved by a locally averaged approximation of the continuity and Navier-Stokes equations.The velocity and pressure fields are obtained by numerically integrating the following set of differential equations: where  represents the interphase momentum transfer per unit volume between the particles and the fluid.The full system is closed with the definition of such term, which in our formulation reads: in which  is the number of particles in the volume ; the forces  and  represent the drag and pressure gradient (or generalized buoyancy) force, respectively; and  is a weight function evaluated at the particle center, which allows the resulting source term in Equation ( 2) to be smoothed out over grid cells.
Our simulations are based on a modelling approach combining the Discrete Element Method for the solid phase and a local average CFD approach for the fluid phase.The

Conventional Approach 2.1. CFD-DEM Modelling for Fluidized Beds
In CFD-DEM, the equations governing the two phases are solved separately, i.e., sequentially.Since the characteristic times can be different, distinct time-steps are used for the CFD and the DEM parts, so that generally the fluid flow remains constant for a number of DEM (smaller) time-steps.
The fluid phase flow is solved by a locally averaged approximation of the continuity and Navier-Stokes equations.The velocity and pressure fields are obtained by numerically integrating the following set of differential equations: where F f p represents the interphase momentum transfer per unit volume between the particles and the fluid.The full system is closed with the definition of such term, which in our formulation reads: in which N p is the number of particles in the volume Ψ; the forces F d and F b represent the drag and pressure gradient (or generalized buoyancy) force, respectively; and w i is a weight function evaluated at the particle center, which allows the resulting source term in Equation ( 2) to be smoothed out over grid cells.
Our simulations are based on a modelling approach combining the Discrete Element Method for the solid phase and a local average CFD approach for the fluid phase.The physical equations governing the motion of the particles and of the fluid are summarized below.
To track the translational and rotational motion of each individual particle in the system, the following equations are solved: where m i , v i , I i and ω i are the i-th particle mass, velocity, moment of inertia and angular velocity, respectively.The summation of external actions includes contact forces, ∑ Nc j=1 F c,ij , the total hydrodynamic force (i.e.drag + pressure gradient), F h,i , gravity, F g,i , and cohesive/adhesive forces F k,i .The cohesive forces allow for the inclusion of different models, such as van der Waals, capillary bridge, and electrostatic effects.
In the rotational direction, the summation is on all torque contributions generated by non-collinear collisions, T c,ij , and the corresponding rolling friction torque, T r,ij , and the fluid-particle torque, T f p,i .

Contact Models
The contact force is computed using the linear spring-dashpot-slider [39] model whose expressions for the normal and tangential component of the force are where the δ's represent the (normal, sub n, and tangential, sub t) displacements between the contacting particles, v their relative velocity components at the contact point, K the spring stiffness constants, η the dashpot damping coefficients and µ the slider friction coefficient.Note that the tangential contribution of the force is capped in magnitude by Coulomb's sliding limit, the rest of the associated energy being dissipated as friction.
The coefficient of restitution, e n , determines the damping coefficient, η n , according to A more adequate representation of the interparticle contact is through Hertz-Mindlin theory [40][41][42], which, neglecting micro-slip on the contacting surfaces (no-slip approximation), is characterized by parameters that can be related to the linear counterpart.
The necessary formulas are reported below: in which the equivalent properties, E eq , R eq , are used in the case of particles of different materials or size [42].The forces can be conveniently expressed as The velocity dependent dissipative terms depend on the parameters η H n and η H t which are related to the corresponding restitution coefficients, e n and e t , respectively, according to: in which Equation ( 10) should be used for K n (see e.g., [43,44]).The same Coulomb limit for the tangential force as in Equation ( 8) applies also to the non-linear contact model.

Hydrodynamic Interaction Models
Models for the drag force or drag coefficient as a function of the slip velocity and voidage are abundant and cover a broad range of slip velocities and voidage values.Available commercial and open-source packages such as Fluent, Star-CCM+, MFiX, openFOAM offer several alternatives.The most common general form is expressed as: in which β contains the dependence on the particle Reynolds number Re p = and voidage ε and can take a complex form.We exemplarily report here the classical model known as Gidaspow [5]: The expressions of other models utilized in simulations of fluidized beds are omitted for brevity, but the interested reader finds useful references in the following (nonexhaustive) list: Di Felice [45], Beetstra or BVK [46], HYS [47,48], Rong et al. [49], Cello et al. [50] and Tang et al. [51].
A set of models has been developed to deal with the presence of polydisperse particles, as DEM requires the force on each individual particle, and semi-empirical expressions can only estimate the average drag force in (portions of) the bed.Their general for is where F d,avg is evaluated with monodisperse models and the specification (or repartition) coefficient γ i distributes such average force across the different particle size classes.It is itself a function of the flow conditions (Re p and ε) and the local polydispersion index in which D avg is Sauter mean diameter , and x k is the volume fraction of the particle size class.Examples of such models are available in [46][47][48][49] and [52].
The pressure gradient, or generalized buoyancy, force is where V p is the particle volume, and ∇p is the gradient of the averaged pressure.
Other contributions such as transient forces (added mass, Basset's history integral), Saffman and Magnus lift, and fluid torque are not commonly employed for fluid beds as the prevalence of dense regions and typical applications to gas-solids systems make them negligible.Moreover, the use of a relatively coarse grid does not make CFD-DEM easily amenable to include them.For reference, some applicable formulations can be found in [53].

Cohesive Forces
Finally, common applications of fluidized beds require the consideration of cohesive and interparticle forces.Examples include fluidization of fine and ultrafine particles, fluid bed granulation and coating, fluidized beds of insulating (e.g., polymer, organic) particles.Detail treatment of the formulation for each case is out of the scope of this review, so a short list of the typical models and recent relevant references is presented below, and additional details will be provided in the context of their coarse graining version:

From Real to Computational Particles (Grains)
Rather than tracking the trajectories of each individual particle in the system, the coarse graining approach offers the attractive idea of lumping together close particles into a computational, representative element.Ideally, the same properties as the original system will be obtained in the coarse-grained system.From a statistical point of view, there are different possibilities to reduce the order of the system.In the context of CG for DEM, the most important ingredient is that the discrete nature of the multi-particle system is maintained.Thus, opposed to considering cell-averaging in a Eulerian grid, the coarsegrained particles are modelled as distinct elements that interact through collision events (among themselves and with walls) and are subjected to the action of fluid drag.
The intended focus of this review is on those methods that clearly rely on a DEM-like approach, for which the equations of motions of the representative and actual particles share the same structure, so that numerical methods and libraries developed for DEM can be readily utilized also for coarse-grained simulations.
The procedure to move from actual to representative particles is the key element of the coarse graining procedure.To keep a consistent nomenclature, we assume here to name as particles the original solid elements, and as grains the representative, or computational coarse-grained elements, corresponding to the parcels (Figure 2).The word "grain" has the advantage that is directly related to the coarse graining procedure it originates from.In addition, a subscript "g" can be used to discriminate grains from particles, for which the subscript "p" will be used.
Two key parameters to define the coarse graining are introduced: the coarse grain factor f and the number of particles per grain n CG , whose definitions are reported Table 1.
Table 1.Definitions of the coarse-graining parameters.

Variable, [Units]
Description coarse grain factor, grain-to-particle size ratio coarse grain number, number of particles in a grain The derivation of relations linking the physical particles' properties to those of the computational grains has followed different approaches.After a short historical context, the following subsections summarize the procedures to link the real particle properties and those of the corresponding coarse grains, separately for the general properties, the contact forces (with both the linear and non-linear Hertz-based approaches), the drag forces and other contributions.Two key parameters to define the coarse graining are introduced: the coarse grain factor  and the number of particles per grain  , whose definitions are reported Table 1.
The derivation of relations linking the physical particles' properties to those of the computational grains has followed different approaches.After a short historical context, the following subsections summarize the procedures to link the real particle properties and those of the corresponding coarse grains, separately for the general properties, the contact forces (with both the linear and non-linear Hertz-based approaches), the drag forces and other contributions.

Context of the Early Coarse-Graining Approaches
In the earliest approaches, the concept of representative particle assumed a statistical representation meaning, mostly in the context of scaling rules and dimensionless groups for fluidized beds.In their attempt, Sakano et al. [18] introduced in DEM the Imaginary sphere model, in which the actual particles were assumed represented by grains 4 to 9.5 times bigger, with appropriately scaled solid density.Collisions were assumed to occur between grains as with the particles, and the drag force was assumed to apply in the same way to the grains as to the particles.Kuwagi et al. [19] introduced the Similar Particle Assembly (SPA) concept, according to which the solids mass and volume were kept constant upon coarse graining, and volume fraction and flow similarity of the particles and the grains was assumed, i.e.,  =  and  =  ( = average particle velocity), respectively.They later used the SPA model to model thermoset particles with a large coarse graining factor,  = 200 (corresponding to 8000 grains for 64 billion particles) [14] and later validated it for bubble size distribution in a 2D gas-fluidized bed of Geldart's group A and D particles [59].Washino et al. [60] also introduced scaling rules considerations by dimensionless groups in DEM simulations of fluid beds, showing attractive computational scaling properties (  = −2.55⋅  ).

Context of the Early Coarse-Graining Approaches
In the earliest approaches, the concept of representative particle assumed a statistical representation meaning, mostly in the context of scaling rules and dimensionless groups for fluidized beds.In their attempt, Sakano et al. [18] introduced in DEM the Imaginary sphere model, in which the actual particles were assumed represented by grains 4 to 9.5 times bigger, with appropriately scaled solid density.Collisions were assumed to occur between grains as with the particles, and the drag force was assumed to apply in the same way to the grains as to the particles.Kuwagi et al. [19] introduced the Similar Particle Assembly (SPA) concept, according to which the solids mass and volume were kept constant upon coarse graining, and volume fraction and flow similarity of the particles and the grains was assumed, i.e., ε g = ε p and v g = v p (v p = average particle velocity), respectively.They later used the SPA model to model thermoset particles with a large coarse graining factor, f = 200 (corresponding to 8000 grains for 64 billion particles) [14] and later validated it for bubble size distribution in a 2D gas-fluidized bed of Geldart's group A and D particles [59].Washino et al. [60] also introduced scaling rules considerations by dimensionless groups in DEM simulations of fluid beds, showing attractive computational scaling properties (log CPUtime = −2.55•logf ).

Compact Grains
In the majority of the coarse graining approaches, the assumptions on the conservation of mass, volume and density of the particles apply in a similar fashion.Grains are assumed to be compact (i.e., non-porous) and account for a given number of particles.The total mass of the solid for particles and grains is set to be the same.The total volume is also typically set to be constant.All grains represent the same number of particles and, for monodisperse systems, are monodisperse.Therefore, the following relations apply: From the relation between the grain and particle size (d g = f d p , see Table 1), it follows that Bearing in mind that the grains are assumed spherical, the constraint of same volume is equivalent to saying that the solids' and gas volume fractions in a given region of the fluidized bed will be automatically maintained, i.e., as the packing degree of spheres is known to be independent of size.

Porous Grains
For fast moving systems, such as risers and turbulent or fast-fluidized systems, the interplay between contact and hydrodynamic forces is shifted towards the latter, or at least dense regions are more limited in space.In addition, the scale of the systems is generally big, so coarse grids for the fluid are also necessary to keep the simulation feasible.For these complex multiphase and highly multiscale dynamical systems, coarse graining of DEM may require special treatments.A specific strategy was introduced by Lu et al. [16], in the context of the Energy-Minimization Multi-Scale (EMMS) method for the simulation of circulating fluidized bed reactors.Coarse grains represent both a given number of actual particles and the void among them and are introduced to represent an intermediate scale between that of the particles, on one side, and that of the typical heterogeneous structures (e.g., clusters), on the other side (Figure 3).The grain voidage is assumed to be the same as the voidage surrounding it.So, for 3D systems, where ε CGP is the void fraction of the grain.Isolated grains are effectively representative if the convective transport mechanism of its elementary particles dominates over the diffusional counterpart.
The scaling correlation for mass is: One of the most detailed treatment and extensive application of the coarse graining procedure with the linear contact model was given by Sakai and coworkers, who introduced it to model a pneumatic conveying line [23] and later applied it for several fluidized bed cases [61][62][63][64][65][66].They based their derivation on imposing the conservation of the typical energy contents of the particle/grain system before and after the collisions.
The grains are assumed to move with the same velocity as the average velocity of the actual particles they represent (Figure 4).More simplistically, grains represent actual particles that move all with the same velocity, which is also equal to the grain velocity, i.e., Processes 2021, 9, x FOR PEER REVIEW 9 of 30

Constant Absolute Overlap Models
One of the most detailed treatment and extensive application of the coarse graining procedure with the linear contact model was given by Sakai and coworkers, who introduced it to model a pneumatic conveying line [23] and later applied it for several fluidized bed cases [61][62][63][64][65][66].They based their derivation on imposing the conservation of the typical energy contents of the particle/grain system before and after the collisions.
The grains are assumed to move with the same velocity as the average velocity of the actual particles they represent (Figure 4).More simplistically, grains represent actual particles that move all with the same velocity, which is also equal to the grain velocity, i.e., All particles are assumed to rotate also with the same angular velocity, which is not the same as the grain's, as shown below.
The total translational and rotational kinetic energy of a grain can be expressed in terms of the sum of the corresponding energies of the represented particles: The equivalence of the rotational kinetic energy requires that the grain rotational velocity is smaller than the value of the actual particles: All particles are assumed to rotate also with the same angular velocity, which is not the same as the grain's, as shown below.
The total translational and rotational kinetic energy of a grain can be expressed in terms of the sum of the corresponding energies of the represented particles: The equivalence of the rotational kinetic energy requires that the grain rotational velocity is smaller than the value of the actual particles: as it can be easily verified by considering that and equating the rotational kinetic energies appearing in Equation (31).Owing to the derivative and integration of the angular velocity, it follows that the same relation applies also to the rotation, θ, and angular acceleration, α: The fact that collisions of grains occur in well-defined instants requires the assumption that all represented particles would collide in that same moment.Consequently, the impact of two grains leads to the consideration of the impact of all the represented particles.Following the results of Equation ( 31), Sakai et al. [23] set the visco-elastic/frictional interaction forces to manifest the same dependence as the energies, yielding a contact force f 3 times greater than that between the actual particles.
Thus, the normal component is: The particle relative overlap and velocity are replaced by the same quantities for the grain, giving: Before discussing the tangential force, let us present the considerations on the integration time-step, set to be for the grains the same as the value for the particles: It is clear from the analogy between the grain stiffness, K ng , and the particle stiffness, K np , that the normal force-displacement dependence in Equation ( 36) is interpreted as a scaling law for the stiffness and damping coefficient of the grains, as follows: The interesting consequence of such parameter scaling is that the coefficient of restitution of the grains is equal to that of the particles, as it can be easily derived from Equation (9).Therefore, also the dissipated energy during collisions is kept constant.
For the tangential component, in the case of adhering contact surfaces: in which the displacement and velocity are to be calculated at the point of contact between the surfaces.Considering that the grain rotation and angular velocity scale inversely and the radius scales directly with f , the displacement and velocity at the contact point of the grain are maintained similar to the ones of the particles.Following the same approach as for the normal force: This can be effectively interpreted as a scaling law for the tangential stiffness and damping coefficient of the grains, as follows: One of the consequences of such choice is that the tangential to normal stiffness ratio remains constant.
Finally, for sliding contact surfaces, Coulomb's limit is scaled simply according to: which corresponds to maintaining the same coefficient of friction between the grain and the particles, i.e., Recently Cai and Zhao [67] noted that for the conservation of the dissipation during sliding the friction coefficient should be scaled according to It is interesting to note that the f 3 scaling of the stiffness leads to the same maximum overlap between colliding grains as for colliding particles, for impacts at the same velocity.Indeed, the oscillation frequency of a linear mass-spring system depends on the ratio between the spring constant and the mass (both scaled with f 3 ), and the maximum overlap is given by the impact velocity divided by this frequency.The elastic energy upon compression is also among the conserved quantities for grains with respect to the particles [62]: A similar approach was adopted by Benyahia and Galvin [68], who assumed a smaller number of special particles (the grains) to be representative of all actual particles.Each grain has a statistical weight, which is equivalent to n CG in the present notation.An illustrative scheme was reported in a later publication (Figure 5).In their analysis on the concept of parcels in the MP-PIC approach (see [22,69]) simulated like DEM grains under shear and riser flow conditions, the authors also found a sensitive dependence on the scaling of the stiffness and damping coefficient, proposing a proportionality with the mass (i.e., f 3 factor) for their scaling.Similarly, a fixed tangential to normal stiffness ratio is used.
Motivated by the differences in the velocity distribution and granular temperature, they addressed the reduced amount of dissipation by introducing a modification to the restitution coefficient of the grains compared to the value for the actual particles.The assumption was that during a grain-grain collision, each particle represented by one grain was thought to collide with all the particles in the other grain.and later further discussed in Lu et al. [70].
The proposed coefficient of restitution for the grains is related to that of the particles by The denominator is defined up to a maximum coarse graining degree n CG for each e p .For the frequent case in which 1, the grain restitution coefficient can be conveniently calculated by which incidentally removes the limitation on the maximum n CG and exhibits an intuitive tendency e g → 0 as n CG → ∞, for any value of the original coefficient of restitution e p .Benyahia and Galvin reported a significant improvement of the results for both the shear flow and riser flow in maintaining consistency upon scaling with n CG [69].The following different scaling of the restitution coefficient, based on the kinetic theory of granular flows, was used in a later work [71], Equation ( 49) produces a slightly slower decrease of e g than Equation (48) at low coarse graining factors f .However, it also sets a maximum degree of coarse graining f < 1   1−e 2 p to yield a real value of e g .
Hilton and Cleary [26,27] adopted a similar approach and coarse graining scaling laws, additionally estimating the computational load to scale with order O f −3 . = 1 +  − 1 .
(49) Equation ( 49) produces a slightly slower decrease of  than Equation (48) at low coarse graining factors .However, it also sets a maximum degree of coarse graining  < to yield a real value of  .
Hilton and Cleary [26,27] adopted a similar approach and coarse graining scaling laws, additionally estimating the computational load to scale with order ( ).

Constant Relative Overlap Models
Starting from a slightly different point, Radl et al. [25] also examined the coarse graining of particles in parcels treated DEM-like for granular jet applications.Their approach differs from the one presented above for two important aspects: the conserved properties and a distinction between the dense and dilute regions, the latter of which requiring a velocity relaxation.The point of departure is the non-dimensional form of the equations of motion of a particle (or grain), for which all quantities are to be conserved when moving from particles to grains.
The main difference concerns the contact parameters, which are shown to scale according to: Such dependence can be explained by the fact that the maximum absolute overlap of the grains is also scaled compared to the original particles.In particular, the relative maximum overlap is kept constant, for example in terms of similar percentage of the radius, i.e., It shall be noted that the elastic energy stored during collisions is preserved also with this scaling.A smaller stiffness compared to the previous scaling affects other dynamic properties of the collision, such as the duration.It turns out that the collision time scales linearly with the coarse graining degree: However, the duration of collisions is generally smaller by orders of magnitude than any other characteristic time, at least for fluidized bed applications.On the other hand,

Constant Relative Overlap Models
Starting from a slightly different point, Radl et al. [25] also examined the coarse graining of particles in parcels treated DEM-like for granular jet applications.Their approach differs from the one presented above for two important aspects: the conserved properties and a distinction between the dense and dilute regions, the latter of which requiring a velocity relaxation.The point of departure is the non-dimensional form of the equations of motion of a particle (or grain), for which all quantities are to be conserved when moving from particles to grains.
The main difference concerns the contact parameters, which are shown to scale according to: Such dependence can be explained by the fact that the maximum absolute overlap of the grains is also scaled compared to the original particles.In particular, the relative maximum overlap is kept constant, for example in terms of similar percentage of the radius, i.e., It shall be noted that the elastic energy stored during collisions is preserved also with this scaling.A smaller stiffness compared to the previous scaling affects other dynamic properties of the collision, such as the duration.It turns out that the collision time scales linearly with the coarse graining degree: However, the duration of collisions is generally smaller by orders of magnitude than any other characteristic time, at least for fluidized bed applications.On the other hand, the integration time step can be increased, with further additional computational savings.This same coarse graining approach was used by Nasato et al. [72].
Similar to Benhyaia and Galvin [68], Radl et al. [25] recognized the need for additional dissipation to take into account intra-grain collisions between the particles, particularly in dilute regions.For this term, the authors derived a modified version of the velocity relaxation model of O'Rourke and Snider [73].Since the formulation requires the definition of several parameters, a complete discussion is omitted, and the original source is recommended to the interested reader.

Contact between Porous Grains
For porous grains to be used in EMMS-DPM [16], contact model during collisions is to be adapted.The soft grains are thought to be able to compress until the elementary particles reach close packing, i.e., ε m f .The grain has an external diameter d g and an internal "hard-core" diameter, defined by in which the voidage is where ε CGP,min depends on EMMS properties such as the minimum cluster diameter and the maximum cluster voidage.
The linear spring-dashpot model is used to model contact between grains in the normal direction, which is detected whenever the hard cores of the grains come into contact (Figure 6).Inter-grain collisions are treated using a slightly simplified version of the linear spring-dashpot model.The normal spring stiffness is set at the minimum value that ensures stability.Normal damping is applied through the coefficient of restitution, which is determined based on the same dissipation of total energy as the real particles and intra-grain collisions, i.e., Processes 2021, 9, x FOR PEER REVIEW 13 of 30 the integration time step can be increased, with further additional computational savings.This same coarse graining approach was used by Nasato et al. [72].Similar to Benhyaia and Galvin [68], Radl et al. [25] recognized the need for additional dissipation to take into account intra-grain collisions between the particles, particularly in dilute regions.For this term, the authors derived a modified version of the velocity relaxation model of O'Rourke and Snider [73].Since the formulation requires the definition of several parameters, a complete discussion is omitted, and the original source is recommended to the interested reader.

Contact Between Porous Grains
For porous grains to be used in EMMS-DPM [16], contact model during collisions is to be adapted.The soft grains are thought to be able to compress until the elementary particles reach close packing, i.e.,  .The grain has an external diameter  and an internal "hard-core" diameter, defined by in which the voidage is where  , depends on EMMS properties such as the minimum cluster diameter and the maximum cluster voidage.The linear spring-dashpot model is used to model contact between grains in the normal direction, which is detected whenever the hard cores of the grains come into contact (Figure 6).Inter-grain collisions are treated using a slightly simplified version of the linear spring-dashpot model.The normal spring stiffness is set at the minimum value that ensures stability.Normal damping is applied through the coefficient of restitution, which is determined based on the same dissipation of total energy as the real particles and intragrain collisions, i.e., The tangential contact is modelled using a purely dissipative force capped to Coulomb's limit for sliding, with a ratio of tangential to normal damping coefficient fixed at √0.8.The time step is set as a fraction of the collision duration of based on the mass of the The tangential contact is modelled using a purely dissipative force capped to Coulomb's limit for sliding, with a ratio of tangential to normal damping coefficient fixed at √ 0.8.The time step is set as a fraction of the collision duration of based on the mass of the grain, i.e., ∆t = 2π m g K ng , benefitting from the computational advantage of this choice.More details are available in the original manuscript.

Summary of Coarse Graining Approaches for the Linear Model
A summary of the coarse graining scaling factors for the linear spring-dashpot-slider model is reported in Table 2.Note that scaling imposing the same coefficient of restitution (i.e., e g = e p ) is equivalent to using the scaled damping coefficient.

Contact Interaction: Hertz-Based Modelling
Non-linear Hertz-based models are common in DEM simulations of quasi-static systems, when interactions at particle-particle and particle-wall contacts are crucial to represent bulk stress-strain relationships.Nonetheless, contacts and dissipation in fluidized beds can become important, for example when cohesive interactions are also involved.As it will be shown, scaling of Hertz-based model is not particularly complicated.
An accurate treatment of coarse graining for Hertz-based contacts was reported by Bierwisch et al. [24] in the framework of a simulation and experimental study of a cavity filling process.It is based on the same postulates regarding the conservation of the density of various energy contributions.A similar gravitational energy density requires the granular phase to have the same total solids mass and volume and each grain has the same density of the particles; this condition also leads to a similar volume fractions for spherical particles.
Imposing the same kinetic energy density requires all the representative particles to move with the same velocity as the grains.The condition to maintain the dissipated energy during collisions constant requires the collision of one grain to dissipate as much as n CG collisions of the represented particles, i.e., If the initial collision velocity is the same, i.e., v 0g = v 0p , and the coefficient of restitution is correspondingly the same, i.e., e g = e p , so that the final velocity is also equal, then the last term of Equation ( 57) is equal to n CG E dp .
The coarse grain model for the contact forces is obtained by examining the equations of motion for the contacting particle.In the present notation, for the normal and tangential direction and neglecting the cohesive forces, they read Processes 2021, 9, 279 15 of 30 and m eq a t = −min µ 2 3 Note that the rotational motion was not considered [24].By dimensional analysis of the equations and the study of relative dimensionless group, conservation of the energy densities leads to the following dependencies: A comparison with the linear model reveals that the non-linear model naturally leads to a constant relative overlap for collisions at the same velocity.It can be observed that maintaining the same physical parameters, i.e., ρ, E, ν, leads to a consistent coarse graining approach.
It should be noted that the formulation of the velocity dependent dissipative force shows some differences with respect to Equation (15).A corresponding dependence of η H n implies a quadratic scaling with the coarse graining factor (see e.g., [44]), i.e., η H n ∝ fˆ2, in order to yield the same constant coefficient of restitution, i.e., e g = e p , similar to previously illustrated linear models.
Overall, the force contributions are shown to scale by a factor of f 2 , a fact that allows also solid stress distribution in dense regions to remain unchanged [24].Indeed, this was also confirmed in the analysis of scaling in simulated experiments of uniaxial compression tests (stress-strain relationship) by Thakur et al. [74].They reported correct scaling of the results with a linear up-scaling of the normal stiffness with the coarse graining degree, which is the case here, as proved by the constant E eq and the linear dependence on f of both R eq and δ n in Equation (10).
In a similar way, Nasato et al. [72] examined the Hertz-based normal contact dependence on particle size and found that it is scale independent.This means that parameters such as Young's modulus and Poisson ratio are set equal.
The consequence of the relative overlap kept constant, i.e., is that the collision duration and integration time-step increase linearly with the coarse graining degree: The same model was used in CFD-DEM simulation of dense medium cyclones by Chu et al. [75].

Hydrodynamic (Drag and Pressure Gradient) Forces
For the drag force, a large consensus points to the need that grains experience the same total force as all the represented particles under the same conditions.For compact coarse grains, this corresponds to assuming the following simple relationship in which β p is to be calculated using the particle properties, not the grain properties, and the fact that v p = v g was used.An implicit assumption of Equation ( 63) is that all represented particles are close enough to one another to experience the same relative velocity and voidage.A similar scaling applies to the pressure gradient force: For the case of non-compact, porous grains, the following modifications were proposed in the context of the EMMS-DPM model: and In Equation (65) the parameter β ph is the conventional coefficient in the drag expression multiplied by a correction factor known as heterogeneity index, β ph = β p H D , defined by in which Re p is the particle Reynolds' number and the constants a, b and c are complex functions of the local voidage ε.A representative plot of the heterogeneity index for a particle of Geldart's group B is shown in Figure 7.
the fact that  =  was used.An implicit assumption of Equation ( 63) is that all represented particles are close enough to one another to experience the same relative velocity and voidage.
A similar scaling applies to the pressure gradient force: For the case of non-compact, porous grains, the following modifications were proposed in the context of the EMMS-DPM model: and In Equation ( 65) the parameter  is the conventional coefficient in the drag expression multiplied by a correction factor known as heterogeneity index,  =   , defined by in which  is the particle Reynolds' number and the constants ,  and  are complex functions of the local voidage .A representative plot of the heterogeneity index for a particle of Geldart's group B is shown in Figure 7.At this point it is worth noting that the use of larger particles (the grains) pushes towards the use of a coarser Eulerian grid for the fluid.In turn, if the cell to grain size ratio is kept constant, increasing the coarse graining degree quickly brings about the issues associated with the use of a large cell size compared to the particle size [77].Therefore, either the need for grid size scaling or sub-grid corrections was recognized and pointed out [78][79][80], as a common requirement also for Eulerian-Eulerian and MP-PIC simulations on coarse grids (Figure 8a).For example, Radl and Sundaresan [78] examined the vertical upflow in periodic domains at different Reynolds' number and solids concentration, later focusing specifically on fluid and particle coarsening for parcel-based simulations [81].Figure 8a shows graphically a typical effect of the use of coarse-grained parcels, which requires fluid grid coarsening.Figure 8b shows the effect of the fluid grid coarsening (filter-to-particle size ratio = 39) on the normalized drag coefficient as a function of porosity and the number of particles per parcel.Remarkably, the dependence on voidage of the At this point it is worth noting that the use of larger particles (the grains) pushes towards the use of a coarser Eulerian grid for the fluid.In turn, if the cell to grain size ratio is kept constant, increasing the coarse graining degree quickly brings about the issues associated with the use of a large cell size compared to the particle size [77].Therefore, either the need for grid size scaling or sub-grid corrections was recognized and pointed out [78][79][80], as a common requirement also for Eulerian-Eulerian and MP-PIC simulations on coarse grids (Figure 8a).For example, Radl and Sundaresan [78] examined the vertical upflow in periodic domains at different Reynolds' number and solids concentration, later focusing specifically on fluid and particle coarsening for parcel-based simulations [81].Figure 8a shows graphically a typical effect of the use of coarse-grained parcels, which requires fluid grid coarsening.Figure 8b shows the effect of the fluid grid coarsening (filter-to-particle size ratio = 39) on the normalized drag coefficient as a function of porosity and the number of particles per parcel.Remarkably, the dependence on voidage of the correction on the drag coefficient (Figure 8b) appears in reasonable agreement with the same dependence of the heterogeneity index in the EMMS model (Figure 7), at least for sufficiently high Reynolds' numbers, although they arise from apparently different points of departure.Both plots highlight the strong need for corrections to the conventional drag models for flow in risers and if significant particle coarse graining is targeted.
correction on the drag coefficient (Figure 8b) appears in reasonable agreement with the same dependence of the heterogeneity index in the EMMS model (Figure 7), at least for sufficiently high Reynolds' numbers, although they arise from apparently different points of departure.Both plots highlight the strong need for corrections to the conventional drag models for flow in risers and if significant particle coarse graining is targeted.

Cohesive Force Modelling
The scaling of cohesive forces is less well established than contact and hydrodynamic forces, so the models appearing will be listed as presented by the authors.On the other hand, scaling for cohesive solids has more severe implications on the accuracy of the computations [82].
In their analysis of the coarse graining using Hertz contacts, Bierwisch et al.
[24] presented also the scaling rule for the Johnson, Kendall, Roberts (JKR) theory of cohesion, considered without hysteresis.The material parameter is the work of adhesion per unit contact area , which is clearly proved to depend on the coarse graining degree leading to a cohesive force scaling with  .A different scaling for the JKR cohesive model was recently introduced by Chen and Elliott [83], who found correct scaling of the cohesive force and surface energy with  and  , respectively.However, the work of adhesion was also scaled to keep the same dissipative restitution (i.e., rebound vs. impact velocity) according to  , which led to the corresponding scaling of the Young's modulus  =  . and the overlap during collisions is obtained to be the same, i.e.,  =  .
In Sakai et al. (2012) [63], the coarse grain scaling on DEM is presented including van der Waals interactions.Similar to the other phenomena, the interparticle cohesion on one grain is applied assuming that all represented particles interact simultaneously with all the particles represented by the close grain.Scaling of the force is initially set to depend on  , in analogy with the other force contributions.Considering that the inter-particle distance, ℎ , may be smaller than the inter-grain distance, ℎ , the following correction is proposed Thus, the overall cohesive force scales with where the Hamaker constant is assumed to be the same for grains as for the particles.The cohesive force does not change the contact force, and the DEM time step is set to be the same as the cohesionless case.

Cohesive Force Modelling
The scaling of cohesive forces is less well established than contact and hydrodynamic forces, so the models appearing will be listed as presented by the authors.On the other hand, scaling for cohesive solids has more severe implications on the accuracy of the computations [82].
In their analysis of the coarse graining using Hertz contacts, Bierwisch et al.
[24] presented also the scaling rule for the Johnson, Kendall, Roberts (JKR) theory of cohesion, considered without hysteresis.The material parameter is the work of adhesion per unit contact area w, which is clearly proved to depend on the coarse graining degree leading to a cohesive force scaling with f 2 .A different scaling for the JKR cohesive model was recently introduced by Chen and Elliott [83], who found correct scaling of the cohesive force and surface energy with f 3 and f 2 , respectively.However, the work of adhesion was also scaled to keep the same dissipative restitution (i.e., rebound vs. impact velocity) according to f 3 , which led to the corresponding scaling of the Young's modulus E g = f 2.5 E p and the overlap during collisions is obtained to be the same, i.e., δ ng = δ np .
In Sakai et al. (2012) [63], the coarse grain scaling on DEM is presented including van der Waals interactions.Similar to the other phenomena, the interparticle cohesion on one grain is applied assuming that all represented particles interact simultaneously with all the particles represented by the close grain.Scaling of the force is initially set to depend on f 3 , in analogy with the other force contributions.Considering that the inter-particle distance, h p , may be smaller than the inter-grain distance, h g , the following correction is proposed Thus, the overall cohesive force scales with f 2 where the Hamaker constant is assumed to be the same for grains as for the particles.The cohesive force does not change the contact force, and the DEM time step is set to be the same as the cohesionless case.In Mokhtar (2012) [59], the Similar Particle Assembly model [20] is integrated with the cohesive liquid bridge force, which is not scaled.Grains have diameter f times that of the particle, with equal density and assuming same velocities.The contact model is LSD.The forces are scaled with the third power of f 3 , but all other properties are kept equal to the original particle case, including the DEM time step.
More recently, Chan and Washino [84] proposed a coarse graining strategy for liquid bridge cohesive interactions between grains for applications in agitated mixers.They based the derivation on the assumption that interfacial interactions such as the liquid bridge cohesion shall scale with f 2 .Therefore, in analogy with the f 2 dependence of the (Hertzbased) contact force, also the liquid bridge force (capillary plus normal and tangential viscous contributions) is calculated by The following additional parameter scaling were shown to be successful.Liquid bridge volume: Separation distance: Rupture distance: The most recent investigation by Tausendschön et al. [85] concerns a detailed analysis of liquid bridge and van der Waals cohesive force scaling in coarse grained CFD-DEM simulations of periodic fluidized systems.In the context of contact scaling based on maintaining constant the relative overlap (see Table 2), starting from three different theoretical bases led essentially to confirm correct scaling if the surface tension and liquid viscosity scale linearly with f and the Hamaker constant scales with f 3 .In addition, they pointed their attention to the significant role that the field smoothing filter plays in hydrodynamics for grain sizes growing similar to or larger than cell sizes.

Computational Savings
As anticipated in the introduction, the main advantages of the coarse graining approach and its most attractive features are the conceptual simplicity, which translates into simple steps to deploy it into a code, and the computational saving compared to classical DEM.For nearly all the approaches presented so far, the changes to the code required to implement coarse graining are very limited.The core of the coarse graining lays in adapting the parameters so that the (same) DEM equations of motion and CFD part provide the results in terms of grains instead of the particles.Therefore, another attractive property of the method is the ease of implementation.
The theoretical scaling of the computational time required to simulate a given time can be estimated assuming an inverse proportionality of the total time on the integration time step (realistic) and on the number of particles (pessimistic).As discussed in Section 3, the key step determining the time step is the contact model: with the constant absolute overlap linear model, the time step with grains is the same as the one with particles; with the constant relative overlap linear model and the Hertz model, the time step increases proportionally with f .The number of particles scales with f 3 .The coarse graining factor represent the ratio between the grain size and the particle size.Therefore, constant absolute overlap and constant relative overlap schemes scale with the 3rd and 4th power of the coarse graining factor, respectively.A plot of these theoretical speed-up trends is shown in Figure 9.In the best case (constant relative overlap), grains twice as big as the particles they represent (n CG = 8) already allow an estimated speed-up = 16.With grains just three times bigger (n CG = 27), the speed-up is greater than 80. Eighty times faster!lected applications from the literature are listed in Table 3.
Considering the one decade since introduction and widespread diffusion, it must be admitted that the question on to what extent the highly attractive savings come at the cost of accuracy is still not fully answered.To investigate validation at extreme coarse graining, a recent application in fluidized beds involved grains with  = 300,000, i.e.,  ≈ 67 [33], but the expected speed-up was not reported.
Figure 9. Scaling performance of coarse graining schemes.The speed-up is the computational gain achievable thanks to the coarse graining to simulate a given time, assuming the computational load scales linearly with the number of particles.Table 3. List of selected applications of coarse graining approaches and corresponding performance improvements.

Reference Investigated CG Range
Typical Performance It is worth mentioning that the above estimate is not optimistic, as a decrease in the number of particles is very likely to produce a higher than linear decrease in CPU time.However, in CFD-DEM simulations part of the time is also spent in the CFD part, with overall savings that can be less striking.
To show examples of performance achievable with coarse graining CFD-DEM, selected applications from the literature are listed in Table 3.
Table 3. List of selected applications of coarse graining approaches and corresponding performance improvements.

Reference
Investigated CG Range Typical Performance Considering the one decade since introduction and widespread diffusion, it must be admitted that the question on to what extent the highly attractive savings come at the cost of accuracy is still not fully answered.To investigate validation at extreme coarse graining, a recent application in fluidized beds involved grains with n CG = 300, 000, i.e., f ≈ 67 [33], but the expected speed-up was not reported.

5.
Applications to Bubbling/Spouted Beds 5.1.Bubbling/Spouted/Liquid Fluidized Beds Initial verification and validation studies have been carried out on small pseudo-2D geometries and bubbling fluidization conditions [61][62][63].More recently, CG-DEM was simulated in lab-scale bubbling beds with immersed tubes [66].Simulations have been reported on the small-scale NETL bubbling bed challenge [88], including also the benefit of hard-sphere model for the collisions, which further considerably improves numerical efficiency.
Coal gasification was studied in a bubbling fluidized bed including heat transfer and heterogeneous chemical reactions, which allowed the influence of the operating parameters to be effectively characterized [89].Using a commercial software, the steam gasification of biomass was studied in a bubbling bed including heat transfer, chemical reaction and particle shrinkage [34,90].The fast pyrolysis of biomass was investigated using a multiscale approach, combining coarse-grain particle scale simulations with reactor scale modelling [91,92]; similarly, coarse grain DEM-CFD combined with reduced-order modelling was used to simulate a pilot-scale updraft coal gasification reactor [93].Bubbling fluidization of a sand-biomass mixture was compared against experiment for degree of mixing and pressure drop (average and fluctuations) [94,95].Different testing methods for measuring solids distribution have been compared and improved by simulating the so-called "travelling fluidized bed" [96].A coarse-graining application to segregation in vibrated fluidized beds can be found in Reference [97].
Very long simulations have been achieved using EMMS-DPM for the methanol-toolefin process [98] (Figure 10), reaching as long as 8 h of simulated time.Similar long runs were obtained with a similar model for the simulation of a continuous compartment fluidized bed, with computation of the residence time distribution of the polydisperse fluidized solids [99].Reprinted from [98], with permission from Elsevier.

Circulating Fluidized Beds and Cyclones
Since the introduction of the coarse grain technique, applications in riser flow and CFBs have started to appear in noticeable number.Several periodic systems were investigated, mostly with the objective to improve understanding, characterize the influence of the coarse graining degree and grid size.Examples of such analyses are [68,78,80,81].
Larger scale riser flow simulations were carried out for energy applications, such as for a chemical looping system [102] and for a fluidized bed carbonator [103,104], in both cases using a commercial software.A complex combination of tools including coarse grain Solids motion in spouted beds was investigated with focus on the solids and velocity distributions, pressure drop and solids mixing [86].The effect of gas inlet in prismatic spouted beds was also investigated successfully [100].As a peculiar application, inter-mittent fluidization during solids discharge of lunar regolith simulants was shown to be successfully simulated [101].
Simulations of a liquid-fluidized reactor, including chemical conversion were reported in [70].

Circulating Fluidized Beds and Cyclones
Since the introduction of the coarse grain technique, applications in riser flow and CFBs have started to appear in noticeable number.Several periodic systems were investigated, mostly with the objective to improve understanding, characterize the influence of the coarse graining degree and grid size.Examples of such analyses are [68,78,80,81].
Larger scale riser flow simulations were carried out for energy applications, such as for a chemical looping system [102] and for a fluidized bed carbonator [103,104], in both cases using a commercial software.A complex combination of tools including coarse grain CFD-DEM was successfully utilized to simulate the multiphase flow of liquid, magnetite fine particles, coal particles, with the formation of an internal air core [105].The combined use of EMMS-DPM has also been extensively used for riser flow, e.g., [16,106].

Physical Models
Since its introduction, the concepts of parcels and grains appeared to be limited to a constant (in time) and uniform (across the system) number of particles per grain.This restriction had severe implications from the point of view of industrial applicability, as real particles are all but monodisperse spherical particles.On the other hand, DEM is notorious for its capability to treat particle individually, so potentially each one different from the rest.In one of their recent efforts, Lu et al. [71] introduced two strategies to deal with coarse graining of a particle system with size distribution.The two strategies differ for the constant quantity of the coarse graining: in strategy SSW (same statistical weight), the grains contain the same number of particles; in the SSP (same size parcel), all the grains share the same size (Figure 11).As a consequence, in the former case, there is a size distribution of the simulated grains; in the latter, there is a number distribution of represented particles by each grain and all grains are monodisperse.The technique was successfully applied to the fluidization of Geldart's group A (FCC) and group B (glass beads) particles and showed that the agreement with pure DEM and experiments is reasonable.Lu et al. [71] conclude that the accuracy of the SSW method at low velocity is higher than that of the SSP method.At higher velocity their accuracy is similar.On the other hand, the monodisperse nature of the grains with the SSP strategy allows for improvements in the computational efficiency, avoiding the limitation on the time step imposed by smaller particles/grains.
Another improvement in coarse-grain simulations concerns accurate procedures for scaling heat transfer parameters.This has been thoroughly investigated by Lu et al. [107], who determined scaling parameters for both the particle-fluid-wall transfer mechanism and the direct contact heat transfer mechanism.The results indicate that: the particle minimum distance and gas layer length must be kept constant in dimensionless terms, i.e., they should scale with the coarse graining degree f ; the dimensionless particle-wall distance scales with f −1 ; so, overall, the particle-fluid-wall transfer mechanism calculated using the grain properties must be scaled by Similarly, for the direct contact mechanism, in which f (e) is an explicit function of the coefficient of restitution.Note that the result for the scaling of the direct contact mechanism assumes the linear spring-dashpot model with the same particle stiffness of the grains as the particles.
Treatment of the erosion of tubes in bubbling fluidized beds was investigated by Zhou and Zhao [108].
Processes 2021, 9, x FOR PEER REVIEW 22 of 30 the SSP method.At higher velocity their accuracy is similar.On the other hand, the monodisperse nature of the grains with the SSP strategy allows for improvements in the computational efficiency, avoiding the limitation on the time step imposed by smaller particles/grains.Another improvement in coarse-grain simulations concerns accurate procedures for scaling heat transfer parameters.This has been thoroughly investigated by Lu et al. [107], who determined scaling parameters for both the particle-fluid-wall transfer mechanism and the direct contact heat transfer mechanism.The results indicate that: the particle minimum distance and gas layer length must be kept constant in dimensionless terms, i.e., they should scale with the coarse graining degree ; the dimensionless particle-wall distance scales with  ; so, overall, the particle-fluid-wall transfer mechanism calculated using the grain properties must be scaled by Similarly, for the direct contact mechanism, in which () is an explicit function of the coefficient of restitution.Note that the result for the scaling of the direct contact mechanism assumes the linear spring-dashpot model with the same particle stiffness of the grains as the particles.
Treatment of the erosion of tubes in bubbling fluidized beds was investigated by Zhou and Zhao [108].

Computational Improvements
Recently proposed computational improvements are aimed either at extending the specific coarse graining technique beyond the previous boundaries or at combining multiple techniques to further speed up the computations, without compromising the accuracy.
Starting from the reduced particle stiffness model [109], the augmented coarse graining method for the CFD-DEM simulations of fluidized beds was proposed by Lin et al. [110], who tested the combined model in bubbling and circulating fluidized beds against experimental averaged solids velocity profiles.The computational benefits of coarse graining ( = 5) is combined with the advantage of the reduced particle stiffness model (reduced stiffness factor = 0.01), reaching a speed-up factor of almost 5000.The actual loss in accuracy requires additional and more extensive quantification.
The multilevel coarse graining technique for DEM introduced by Queteshiner et al. [44,111] aims at making it possible to adapt the coarse graining level to different regions, so that a coarser realization is used in regions whose behavior is simpler to represent and

Computational Improvements
Recently proposed computational improvements are aimed either at extending the specific coarse graining technique beyond the previous boundaries or at combining multiple techniques to further speed up the computations, without compromising the accuracy.
Starting from the reduced particle stiffness model [109], the augmented coarse graining method for the CFD-DEM simulations of fluidized beds was proposed by Lin et al. [110], who tested the combined model in bubbling and circulating fluidized beds against experimental averaged solids velocity profiles.The computational benefits of coarse graining ( f = 5) is combined with the advantage of the reduced particle stiffness model (reduced stiffness factor = 0.01), reaching a speed-up factor of almost 5000.The actual loss in accuracy requires additional and more extensive quantification.
The multilevel coarse graining technique for DEM introduced by Queteshiner et al. [44,111] aims at making it possible to adapt the coarse graining level to different regions, so that a coarser realization is used in regions whose behavior is simpler to represent and a finer level is used in spatially confined regions of interest.Different levels are coupled by exchanging volume-average flow properties.One application is the discharge of solids from a bin, taken as reference test for comparison.It is a method with promising potential for improved accuracy and likely to be useful also in fluidized bed applications.On the other hand, the improved accuracy and representation capability comes at the cost of the computational efficiency, which deteriorates compared to pure coarse graining and further work is needed to recover as much as possible of this attractive property.
A very significant improvement in speed for computations is obtained by combining coarse graining and hard-sphere event-driven models for collisions, without compromising the accuracy of the predictions in fluidized beds, according to Lu et al. [87,88] (Figure 12).Together with coarse-graining, a corrected hard-sphere model specifically derived by the authors, formed the CGHS (coarse grain hard-sphere) model, capable of providing an overall impressive speed-up of 24 using f = 2.The solver was fast enough to require computational times comparable with those of MP-PIC, known to be less accurate but much faster.The actual comparison with MP-PIC simulations on the same system (smallscale fluidized bed) revealed a clear superiority of the discrete CGHS model.Indeed, the conclusions of the authors are that the CGHS combinations provide an accurate and inexpensive way to simulate fluidized beds.
Finally, other dynamic coarse graining strategies were proposed by Wang et al. [112] with particle-group crushing and polymerization.The method requires additional corroboration and more extensive validation.
12).Together with coarse-graining, a corrected hard-sphere model specifically derived by the authors, formed the CGHS (coarse grain hard-sphere) model, capable of providing an overall impressive speed-up of 24 using  = 2.The solver was fast enough to require computational times comparable with those of MP-PIC, known to be less accurate but much faster.The actual comparison with MP-PIC simulations on the same system (smallscale fluidized bed) revealed a clear superiority of the discrete CGHS model.Indeed, the conclusions of the authors are that the CGHS combinations provide an accurate and inexpensive way to simulate fluidized beds.
Finally, other dynamic coarse graining strategies were proposed by Wang et al. [112] with particle-group crushing and polymerization.The method requires additional corroboration and more extensive validation.

Open Problems
To really reach sufficient maturity and industrial relevance like it is now for CFD, the attractive computational gain achieved with coarse-grain CFD-DEM still requires substantial effort in terms of further development and validation.Not only are simulations to be confronted with full CFD-DEM results as in many of the articles cited here but also directly with experimental data, as the number of modelling "layers" and corresponding parameters is high.The coarse graining procedure enhances the importance of accurate parameters and requires closed-form procedures for particle-to-grain property scaling.Therefore, it would be highly desirable to compare simulated results with sufficiently simple experiments, whose measured data relate as directly as possible to the microscopic mechanisms that originate them.
To conclude, the following topics, insufficiently covered in previous studies, are proposed as future objectives of interest:

•
Coarse graining of non-spherical particles; • Rotational motion and rolling friction for applications e.g., in cyclones; Copyright (2017) American Chemical Society.

Open Problems
To really reach sufficient maturity and industrial relevance like it is now for CFD, the attractive computational gain achieved with coarse-grain CFD-DEM still requires substantial effort in terms of further development and validation.Not only are simulations to be confronted with full CFD-DEM results as in many of the articles cited here but also directly with experimental data, as the number of modelling "layers" and corresponding parameters is high.The coarse graining procedure enhances the importance of accurate parameters and requires closed-form procedures for particle-to-grain property scaling.Therefore, it would be highly desirable to compare simulated results with sufficiently simple experiments, whose measured data relate as directly as possible to the microscopic mechanisms that originate them.
To conclude, the following topics, insufficiently covered in previous studies, are proposed as future objectives of interest:

Conclusions
Starting from the beginning of the 2000s, a new class of methods to expand the applicability of CFD-DEM have established in the field for the simulation of fluidized beds.Thanks to their significant computational advantages, in the last decade, the coarsegraining approach has grown to become a relatively widespread technique.The present review summarizes the essential ingredients of the coarse graining method, covering in detail: contact forces (linear and Hertz-based), hydrodynamic interactions and cohesive forces for both compact and porous grains.More than one possible approach is found, depending whether the absolute or relative displacement is kept constant, and discussed.As a byproduct of this process, new links between the scaling with the linear and Hertzbased contact model are highlighted.Scaling of the most important hydrodynamic forces, drag and pressure gradient, in homogeneous systems is obtained relatively easily and ensures conservation of the pressure drop for uniform systems (e.g., the dense phase).Where heterogenous structures are dominating (e.g., clustering in risers), the use of a coarse grid for the CFD part often requires additional corrections and adjustments on the drag law.Scaling of interparticle cohesive forces depends on the force model and can be different in the case of Hamaker's model, JKR theory or capillary interactions due to liquid bridges.
The advantages in terms of computational savings is quantified, depending on the coarse graining choice, i.e., constant absolute overlap vs. constant relative overlap models, showing 3rd or 4th power scaling, respectively.Some applications in the field of bubbling and circulating fluidized beds are illustrated, to provide a picture of the state-of-the art and showcase the explored possibilities.New extensions for accuracy improvement (e.g., handling polydisperse solids, scaling of heat transfer and adaptive coarse graining) and further computational saving (combined coarse graining hard-sphere, CGHS) are discussed.
In terms of model fundamentals, this review shows that the basic equations are ready for implementation for most applications.Thus, CGDEM simulations of fluidized bed systems can benefit from unparallel speedup factors.The loss of accuracy is generally limited, for example, if heterogeneities are not dominating, and there are no key particlewall interactions, such as in highly complex geometries.In these cases, the actual loss is yet to be fully quantified, particularly if the coarse graining degree is high.Therefore, expected progress is in the area of validation, possibly with directly comparable experiments (see e.g., [35]) and attempts to quantify the uncertainty associated with coarse graining at different scales (see e.g., [113]).Moreover, developments to include new, locally adaptive schemes would be particularly useful.Overall, the increasing impact of the coarse grained CFD-DEM in recent years is expected to grow further and lead to widespread adoption, both in fluidized bed research and industrial applications.

Conflicts of Interest:
The authors declare no conflict of interest.The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript; or in the decision to publish the results.

Figure 1 .
Figure 1.Number of indexed articles appeared in the last 15 years related to coarse-graining and DEM.Source: data extracted by the database Scopus ® with the following research key: TITLE-ABS-KEY (("coarse-grain*") and ("DEM" or "Discrete Element M*")).

Figure 1 .
Figure 1.Number of indexed articles appeared in the last 15 years related to coarse-graining and DEM.Source: data extracted by the database Scopus ® with the following research key: TITLE-ABS-KEY (("coarse-grain*") and ("DEM" or "Discrete Element M*")).

Figure 2 .
Figure 2. Elements tracked in DEM simulations: (a) elementary particles in a given volume; (b) coarse-grain particles ("grains") in the same volume.The level of coarse-graining and the characteristics of the link between representative and original particles depend on the specific approach.

Figure 2 .
Figure 2. Elements tracked in DEM simulations: (a) elementary particles in a given volume; (b) coarsegrain particles ("grains") in the same volume.The level of coarse-graining and the characteristics of the link between representative and original particles depend on the specific approach.

Figure 3 .
Figure 3. Definition of porous coarse grains in the Energy-Minimization Multi-Scale (EMMS)-DPM model.Reprinted from Reference [16], with permission from Elsevier.

Figure 3 .
Figure 3. Definition of porous coarse grains in the Energy-Minimization Multi-Scale (EMMS)-DPM model.Reprinted from Reference [16], with permission from Elsevier.

Figure 3 .
Figure 3. Definition of porous coarse grains in the Energy-Minimization Multi-Scale (EMMS)-DPM model.Reprinted from Reference [16], with permission from Elsevier.

Figure 4 .
Figure 4. Actual particles and representative grains (CGM) according to the approach of Sakai and Koshizuka [23] for (a) translational and (b) rotational motion.Adapted from Reference [23], with permission from Elsevier.

Figure 4 .
Figure 4. Actual particles and representative grains (CGM) according to the approach of Sakai and Koshizuka [23] for (a) translational and (b) rotational motion.Adapted from Reference [23], with permission from Elsevier.

Figure 5 .
Figure5.The coarse graining concept adopted by Lu et al.[70] for the simulation of liquid fluidized bed reactor.The statistic weight  is equivalent to the number of particles per grain  .Reprinted with permission from Reference[70].Copyright (2016) American Chemical Society.

Figure 5 .
Figure 5.The coarse graining concept adopted by Lu et al. [70] for the simulation of liquid fluidized bed reactor.The statistic weight W is equivalent to the number of particles per grain n CG .Reprinted with permission from Reference [70].Copyright (2016) American Chemical Society.

Figure 6 .
Figure 6.Contact model between grains in the EMMS-DPM approach: (a) actual particles; (b) representative grains in contact and linear spring-dashpot model.Note that  =  in the present notation.Reprinted from [16], with permission from Elsevier.

Figure 6 .
Figure 6.Contact model between grains in the EMMS-DPM approach: (a) actual particles; (b) representative grains in contact and linear spring-dashpot model.Note that d CGP = d g in the present notation.Reprinted from [16], with permission from Elsevier.

Figure 7 .
Figure 7. Dependence of the heterogeneity index  of the EMMS model on the Reyonlds' number and voidage for a particle of Geldart's group B. Reprinted from [76], with permission from Elsevier.

Figure 7 .
Figure 7. Dependence of the heterogeneity index H D of the EMMS model on the Reyonlds' number and voidage for a particle of Geldart's group B. Reprinted from [76], with permission from Elsevier.

Figure 8 .
Figure 8.Effect of the fluid and particle coarsening on the drag coefficient at different coarse graining degree.(a) Graphical example of the fluid-particle coarsening sequence leading to large cell sizes compared to the actual particles; (b) effect of the solid volume fraction ( = 1 − ) and number of particles per grain ( =  in the present notation) on the normalized filtered drag coefficient .Reprinted from [81], with permission from Elsevier.

Figure 8 .
Figure 8.Effect of the fluid and particle coarsening on the drag coefficient at different coarse graining degree.(a) Graphical example of the fluid-particle coarsening sequence leading to large cell sizes compared to the actual particles; (b) effect of the solid volume fraction (φ = 1 − ε) and number of particles per grain (n p = n CG in the present notation) on the normalized filtered drag coefficient β.Reprinted from [81], with permission from Elsevier.

[ 19 ] 29 Figure 9 .
Figure 9. Scaling performance of coarse graining schemes.The speed-up is the computational gain achievable thanks to the coarse graining to simulate a given time, assuming the computational load scales linearly with the number of particles.

Figure 10 .
Figure 10.Simulation of the methanol-to-olefin (MTO) process in a bubbling fluidized bed for 8 h of simulated time, with details on the solids flow (a), age (b) and methanol reaction rate (c).Reprinted from [98], with permission from Elsevier.

Figure 10 .
Figure 10.Simulation of the methanol-to-olefin (MTO) process in a bubbling fluidized bed for 8 h of simulated time, with details on the solids flow (a), age (b) and methanol reaction rate (c).Reprinted from [98], with permission from Elsevier.

Figure 12 .
Figure 12.Schematic concept of the coarse-grain hard-sphere (CGHS) model and performance improvement comparison with other approaches for the simulation of a small-scale fluidized bed.Reprinted with permission from Reference [88].Copyright (2017) American Chemical Society.

Figure 12 .
Figure 12.Schematic concept of the coarse-grain hard-sphere (CGHS) model and performance improvement comparison with other approaches for the simulation of a small-scale fluidized bed.Reprinted with permission from Reference [88].Copyright (2017) American Chemical Society.

Author
Contributions: A.D.R., conceptualization, methodology, resources, writing-review and editing, formal analysis, supervision; E.S.N., investigation, data curation, writing-original draft preparation; F.P.D.M., supervision, validation, writing-review and editing.All authors have read and agreed to the published version of the manuscript.Funding: This research was partially funded on a research contract by Basell Poliolefine Italia (IT), whose support is gratefully acknowledged.Institutional Review Board Statement: Not applicable.Informed Consent Statement: Not applicable.Data Availability Statement: Not applicable.

Table 1 .
Definitions of the coarse-graining parameters.

•
Rotational motion and rolling friction for applications e.g., in cyclones; • Further testing on multi-component polydisperse particles, with coarse-graining of polydisperse drag; • Hydrodynamic interactions in dilute regions (collisional regime) at extreme coarse graining degrees and with coarse grids; • Further development of adaptive coarse graining; • Scaling for mass transfer, chemical reactions and other cohesive interactions, e.g., triboelectric charging.