Coupled CFD-DEM Model Overview

• Coupled CFD-DEM models are Eulerian-Lagrangian frameworks that simulate fluid flow as a continuum while explicitly tracking individual particles.
• They use unresolved and resolved approaches to map particle properties onto CFD grids, ensuring accurate fluid-particle interaction via two-way coupling.
• Advanced implementations attach thermodynamic and reactive state variables to particles, enabling detailed modeling of heat, mass transfer, and phase change.

A coupled CFD-DEM model is an Eulerian-Lagrangian multiphysics framework for particulate flows in which the fluid phase is solved as a continuum and the particulate phase is tracked explicitly as individual entities. In the classical “combined continuum and discrete model” sense, the fluid phase is treated with CFD on an Eulerian mesh, while the solid phase is represented with DEM; the coupling is realized through exchange terms for mass, momentum, energy, and, in reactive systems, species. In the XDEM formulation, this architecture is extended by attaching thermodynamic state variables to each particle, so that particles are treated not only as moving bodies but also as reactive, heat- and mass-transferring micro-systems (Peters et al., 2018).

1. Conceptual basis and model classes

The defining feature of a coupled CFD-DEM model is the simultaneous resolution of continuum-scale fluid transport and particle-scale mechanics. The fluid is represented by CFD or related Eulerian formulations, whereas the solid phase is represented by DEM, with each particle advanced in a fully Lagrangian manner. This arrangement preserves particle identity, morphology, and state, in contrast to Eulerian-Eulerian methods, which are efficient but lose particle-scale information through averaging (Peters et al., 2018).

Within this general framework, two broad realizations recur in the literature. In unresolved, volume-averaged CFD-DEM, the CFD mesh is much larger than the particles, and particle properties are mapped to each fluid cell by averaging over the particles inside it. The dry powder inhaler study is explicit that this is an unresolved, volume-averaged CFD-DEM coupling (Ponzini et al., 27 Aug 2025). In resolved CFD-DEM, the local fluid momentum balance is evaluated over the set of grid cells intersecting each particle, so that fluid-particle forces are obtained directly from resolved pressure and viscous fields rather than from empirical drag laws alone (Fonceca et al., 2021).

The practical distinction between these classes is methodological rather than conceptual. Both retain two-way coupling. In both, the fluid pushes on the solids through drag, pressure-related forces, buoyancy or hydrostatic pressure, viscous interaction, and, depending on the application, heat and mass transfer. The particles, in turn, modify the fluid through porosity, momentum-exchange source terms, and, in thermal or reactive problems, through heat release, heat absorption, evaporation, devolatilization, or species production (Peters et al., 2018).

2. Governing equations and interphase exchange

On the CFD side, one common formulation is incompressible Newtonian flow with gravity as a body force,

$$\rho_f \frac{d\mathbf{U}}{dt} + \rho_f(\mathbf{U}\cdot\nabla)\mathbf{U} = -\nabla p + \mu \nabla^2 \mathbf{U} + \rho_f \mathbf{g},$$

subject to

$\nabla \cdot \mathbf{U} = 0.$

This resolved formulation is advanced with a PISO algorithm in the cited implementation (Fonceca et al., 2021). In more general Eulerian volume-averaged multi-fluid formulations, the important point is that the coupling to the particulate phase enters through the source/sink term $S(\phi)$. Positive values correspond to transfer from solid to fluid and vice versa; $m'$ couples mass exchange, $F$ couples momentum exchange, $Q$ couples heat transfer, and reactive cases add species source terms (Peters et al., 2018).

On the DEM side, particle motion follows Newton-Euler dynamics. A representative formulation is

$$m_i \ddot{\mathbf{x}}_i = m_i \mathbf{g} + \mathbf{F}_{i,f} + \sum_{j=1}^{N_p} \mathbf{F}_{i,j} + \sum_{w=1}^{N_w} \mathbf{F}_{i,w},$$

and

$$I_i \frac{d\boldsymbol{\omega}_i}{dt} = \sum_{j=1}^{N_p} \mathbf{r}_{i,j}\times \mathbf{F}^{\perp}_{i,j} + \mathbf{T}_{i,f},$$

where the right-hand side includes fluid-particle interaction, particle-particle contacts, and particle-wall contacts (Fonceca et al., 2021). In soft-sphere DEM, contact forces, gravity, external forces, and torques are resolved per particle (Peters et al., 2018).

In resolved coupling, the particle-fluid force can be computed cellwise over the set $\mathcal{T}_i$ of grid cells intersecting particle $i$:

$\nabla \cdot \mathbf{U} = 0.$0

The factor $\nabla \cdot \mathbf{U} = 0.$1 multiplies the pressure-gradient term and accounts for the particle shape or solid occupancy inside the cell. In the cited method, gravity is applied in the fluid momentum equation, so the pressure field contains the hydrostatic contribution and the particle does not need an explicit buoyancy correction of the form $\nabla \cdot \mathbf{U} = 0.$2; buoyancy emerges from the pressure field (Fonceca et al., 2021).

The coupled exchange is fundamentally two-way. Fluid-to-particle coupling includes drag, buoyancy, pressure, convective heat and mass transfer, radiation, and turbulence effects. Particle-to-fluid coupling includes momentum sinks or sources, mass transfer from evaporation or reaction, and heat release or absorption (Peters et al., 2018).

3. Numerical realizations, mapping strategies, and discretization choices

Several numerical architectures coexist under the CFD-DEM designation.

Coupling realization	Characteristic	Representative use
Unresolved volume-averaged	CFD cells larger than particles; cellwise averaging of particle data	DPI aerosolization and many packed or fluidized systems
Resolved	Local force evaluation on cells intersecting each particle	Single-particle and confined settling benchmarks
Dual-grid multiscale	Coarse grid for interaction, fine grid for fluid solve	Three-phase flows and high-gradient gas-solid systems

In mono-scale unresolved CFD-DEM, the same grid is often used both to average Lagrangian particle data to Eulerian fields and to solve the fluid equations. Dual-grid multiscale methods separate these roles: a coarse or bulk grid is used for Eulerian projection and fluid-particle interaction, while a fine grid is used for the fluid solution. This separation was introduced to improve grid convergence and accuracy when the grid appropriate for fluid resolution is not the same as the grid appropriate for DEM coupling (Pozzetti et al., 2018).

Porosity mapping is central because local porosity controls drag, momentum exchange, and phase distribution. In packed-bed examples, porosity may be computed per CFD cell using a cell-based algorithm (Peters et al., 2018). However, direct grid-based mapping is strongly grid dependent in unresolved CFD-DEM. A two-grid Voronoi-based coarse-graining method therefore computes a particle-scale local porosity using a radical Voronoi tessellation and maps those values to the CFD mesh through a fine regular point cloud. The reported outcome is a smooth porosity field, reasonable grid independence, and improved particle-scale drag prediction for polydisperse systems (Che et al., 2022).

A more recent semi-resolved strategy inserts an intermediate coarse-graining layer between DEM particles and CFD grids. Particles are first converted into smooth coarse-grained continuum fields via a multi-layer Fibonacci point cloud, independent of the fluid grids; those fields are then projected onto the fluid grid and interpolated back to particles. The purpose is to reduce grid-resolution dependence, suppress oscillations when particles move across grid boundaries, and capture pore-pressure feedback in dense granular systems (Liu et al., 11 Jun 2025).

Discretization choices also vary. A monolithic finite-element CFD-DEM solver in Lethe supports dynamically load-balanced parallelization and high-order schemes, including BDF2, BDF3, and Q2-Q1 elements, while using the same mesh for fluid and particles. A plausible implication is that higher-order FEM makes unresolved CFD-DEM more practical when local volume-averaging requires cells larger than particles (Geitani et al., 2022). Mesh-based FVM and mesh-less SPH have also both been coupled with DEM under locally averaged Navier-Stokes equations; in the comparative study, differences between DEM-SPH and DEM-FVM were found mostly due to differences in computed local fluid fractions (Markauskas et al., 2016).

4. Multiphysics extensions: heat transfer, reaction, phase change, and deformation

A major extension beyond classical DEM is the attachment of internal thermodynamic states to each particle. In XDEM, this includes temperature distributions inside the particle, species distributions in porous particles, reaction progress or conversion, phase change such as melting or drying, and spatially resolved transport coefficients and properties. The review emphasizes that the particle is no longer a rigid object with only position and velocity, but a resolved micro-reactor or micro-heat-exchanger (Peters et al., 2018).

For porous or reactive particles, XDEM uses one-dimensional transient balances inside the particle. A representative particle energy equation is

$\nabla \cdot \mathbf{U} = 0.$3

with $\nabla \cdot \mathbf{U} = 0.$4 containing reaction enthalpy terms or, in melting, latent-heat terms (Peters et al., 2018). At the surface, the particle exchanges heat and species with the surrounding CFD phase through convective, radiative, and conductive terms; radiative and conductive contributions may also be written as inter-particle sums, so the framework includes particle-particle heat exchange and contact effects (Peters et al., 2018).

High-temperature applications sharpen the distinction between mechanical and thermodynamic coupling. In the experimental blast furnace study, XDEM is used to model complex gas-solid flow during the softening process of pre-reduced iron ore pellets. The model includes particle deformation, displacement, temperature, and gas pressure drop and flow; particle overlap increases as Young’s modulus decreases with temperature, reducing local void space and changing gas pathways. To capture the high gas velocity inlet accurately, a dual-grid multi-scale approach is applied (Baniasadi et al., 2018).

Recent thermal models also modify how conduction is represented. An advanced heat transfer model for Eulerian-Lagrangian gas-solid simulations retains Lagrangian DEM for particle motion but treats heat transfer in the solid phase in an Eulerian framework using local void fraction rather than explicit particle contacts. In that model, solid temperature is not computed from particle-by-particle contact conduction and is not influenced by particle diameter or stiffness; the proposed formulation is explicitly presented as independent of spring stiffness and compatible with coarse-grained DEM (Imatani et al., 16 Jun 2025).

5. Engineering domains, validation practice, and representative applications

The documented application space is broad. The XDEM review presents CFD-DEM and XDEM as enabling packed-bed reactors and trickle beds, blast furnace and ironmaking systems, thermal conversion of biomass and solid fuels, melting and phase change, powder metallurgy, fluidized beds and granular reactors, and turbulent three-phase flows (Peters et al., 2018). The recurring advantage across these domains is that local interactions—particle size distribution, porosity heterogeneity, contact networks, reaction fronts, and evolving thermodynamic state—are retained rather than averaged away (Peters et al., 2018).

Sediment transport is a canonical example of high-fidelity CFD-DEM usage. In SediFoam, current-induced transport is modeled across flat bed in motion, small dune, vortex dune, and suspended transport regimes. A central claim of that study is that explicitly retaining particle volume fraction in the fluid equations and drag law changes dense near-bed flow and sediment transport rate in an important way, leading to larger predicted transport rates and better agreement with experiments than prior CFD-DEM sediment studies that did not include those effects in the same way (Sun et al., 2015).

Dry powder inhalers illustrate a different validation culture, combining fluid diagnostics and powder diagnostics. For the NextHaler configuration, the coupled CFD-DEM model is validated against previously published optical experiments. The fluid phase is compared to outlet velocity profiles, 2D PIV velocity maps, and pressure drop profiles; the powder phase is compared to fast-camera plume imaging and optical obscuration signals. The authors report quantitatively good agreement for the time evolution and spatial structure of the plume, and the calibrated DEM model reproduces the carrier emission curve best when the restitution coefficient is around $\nabla \cdot \mathbf{U} = 0.$5 and both sliding and rolling friction coefficients are in the range $\nabla \cdot \mathbf{U} = 0.$6 to $\nabla \cdot \mathbf{U} = 0.$7 (Ponzini et al., 27 Aug 2025).

Fracture, erosion, and deformable porous media extend the framework from particle transport to particle-supported structure. In erosive wear, a conventional solid object is replaced by bonded spherical particles, and fluid, free particles, and the bonded wall interact in a three-way coupled manner; the cited model reproduces classical W-shape and U-shape erosion patterns while allowing surface geometry to evolve through bond failure and particle detachment (Nguyen et al., 15 Feb 2025). For non-Newtonian polymer injection in granular media, the coupled CFD-DEM model reproduces experimental trends in peak pressure and pressure-to-stress ratio, and reports fracture initiation above $\nabla \cdot \mathbf{U} = 0.$8 and $\nabla \cdot \mathbf{U} = 0.$9 (Kazidenov et al., 2024). In a related but distinct DEM-based pore-network framework for deformable porous media, fully coupled hydro-mechanical multiphase flow is validated against Hele-Shaw experiments and used to extract viscous and capillary pressures, energies, contact forces, and flow resistances at grain scale (Dai et al., 2024).

6. Computational scaling, assumptions, misconceptions, and limitations

Computational cost remains one of the main drawbacks of CFD-DEM coupling, especially when particle-scale resolution, thermal physics, or large multiphase systems are included (Peters et al., 2018). This has motivated substantial work on coupling-aware parallelization. A co-located partition strategy places overlapping CFD and DEM subdomains on the same MPI process whenever possible, reducing inter-physics data exchange to direct memory access; in the reported large benchmark, coupled execution on 280 cores used less than 0.1% of the time for inter-physics data exchange (Pozzetti et al., 2018). A parallel dual-grid multiscale strategy further separates DEM and coarse-grid co-location from independent partitioning of the fine CFD grid, and maintains good parallel performance when operated over 1000 processes (Pozzetti et al., 2018).

Time integration of the granular phase is another bottleneck because collision timescales vary strongly across the domain. An adaptive ORB-based local timestepping method therefore advances only collision-prone subsets of particles with very small timesteps while dilute regions use larger timesteps. In the reported mixed dense-dilute cases, the method is 2–3X faster than traditional explicit methods while maintaining the same level of accuracy (Sitaraman et al., 2018).

Model assumptions are equally important. Common simplifications include spherical particles, mono-sized systems, one-dimensional radial internal transport models, preserved sphericity during melting, and empirical or semi-empirical transfer laws (Peters et al., 2018). Some recent work is explicitly aimed at these weaknesses. It is therefore a misconception that a finer CFD grid alone guarantees a better unresolved CFD-DEM simulation: direct grid-based mapping can remain strongly grid dependent, can become unstable as particles move across fluid grids, and can fail to capture pore fluid pressure effects in very dense granular systems (Liu et al., 11 Jun 2025). It is also a misconception that unresolved CFD-DEM is purely a low-order methodology: monolithic high-order implementations now exist, although local averaging still imposes restrictions on particle-to-cell size ratios and Model B is reported to be stiffer and harder to converge than Model A in the cited high-order study (Geitani et al., 2022).

Taken together, these results suggest a stable encyclopedic characterization of the coupled CFD-DEM model. It is not a single algorithm but a family of Eulerian-Lagrangian formulations whose common core is bidirectional exchange between a continuum fluid solver and a particle-resolved mechanical solver. The principal axis of development has been from classical momentum coupling toward richer treatments of porosity mapping, thermodynamics, reaction, deformation, multiscale discretization, and scalable parallel execution, while the main unresolved issues remain empirical closures, computational expense, and the treatment of dense, highly polydisperse, turbulent, or strongly non-spherical systems (Peters et al., 2018).