Eulerian Granular Model
Theory and Physics
Overview
Professor, what is an Eulerian granular model? Does it treat particles as a continuum?
Exactly. The Eulerian Granular Model treats a collection of powders or granules as a pseudo-continuum called the "granular phase" and solves it together with the gas phase within the Euler-Euler framework. It is widely used for gas-solid two-phase flows such as fluidized beds, pneumatic conveying, cyclones, and powder mixing.
How is it different from DEM-CFD?
DEM-CFD tracks individual particles discretely, but for industrial scales involving millions to billions of particles, the computational cost becomes enormous. The Eulerian granular model treats particle groups as a continuum, allowing solutions for large-scale systems within realistic computation times. However, it cannot directly handle individual particle contact forces or shape effects.
Governing Equations
What are the equations for the granular phase?
To describe the motion of the granular phase, we use KTGF (Kinetic Theory of Granular Flow). The continuity equation for the solid phase is the same as the usual Euler-Euler method, but the characteristic feature is using constitutive relations derived from KTGF for the solid phase stress tensor.
Is KTGF the particle version of the kinetic theory of gases?
Exactly. It characterizes particle velocity fluctuations with the Granular Temperature $\Theta_s$.
The transport equation for Granular Temperature is as follows.
The right-hand side terms are, in order: generation by shear, diffusion ($\kappa_s$ is the diffusion coefficient), dissipation due to inelastic collisions ($\gamma_s$), and gas-solid energy exchange ($\phi_{gs}$).
Solid Phase Constitutive Relations
How are solid phase pressure and viscosity determined?
The solid phase pressure $p_s$ is obtained from Granular Temperature and volume fraction. The model by Lun et al. (1984) is representative.
Here, $e_{ss}$ is the coefficient of restitution between particles, and $g_0$ is the radial distribution function. $g_0$ increases sharply as the volume fraction approaches the maximum packing fraction $\alpha_{s,max}$ (about 0.63), representing contact pressure in densely packed states.
| Constitutive Relation | Model Example | Physical Quantity |
|---|---|---|
| Solid Phase Pressure | Lun, Syamlal-O'Brien | $p_s(\alpha_s, \Theta_s)$ |
| Solid Phase Viscosity | Gidaspow, Syamlal | $\mu_s(\alpha_s, \Theta_s)$ |
| Solid Phase Bulk Viscosity | Lun et al. | $\lambda_s$ |
| Frictional Stress | Schaeffer, Johnson-Jackson | Stress in dense packing region |
The Physics of an Hourglass—Granular Matter is Neither "Fluid" Nor "Solid"
Why does sand in an hourglass appear to "flow," yet remains stationary on a sandpile slope? The Eulerian granular model holds the key to answering this question. Granular matter is a "third state" requiring both hydrodynamic continuum equations and solid-mechanical elastic pressure. The concept of granular temperature is an ingenious idea that treats particle velocity fluctuations as temperature, analogous to molecular thermal motion, derived from the kinetic theory of gases by Jenkins & Richman in the 1980s. Without this theory, it would be difficult to predict even the particle circulation velocity in a CFB boiler within an order of magnitude.
Physical Meaning of Each Term
- Temporal Term $\partial(\rho\phi)/\partial t$: Imagine the moment you turn on a faucet. At first, water comes out spluttering and unstable, but after a while, it becomes a steady flow, right? This "period of change" is described by the temporal term. The pulsation of blood flow from a heartbeat, or the flow fluctuation each time an engine valve opens and closes—all are unsteady phenomena. So what is steady-state analysis? It looks only at "after sufficient time has passed and the flow has settled down"—meaning setting this term to zero. Since computational cost drops significantly, trying a steady-state solution first is a basic CFD strategy.
- Convection Term $\nabla \cdot (\rho \mathbf{u} \phi)$: What happens if you drop a leaf into a river? It gets carried downstream by the flow, right? This is "convection"—the effect where fluid motion transports things. Warm air from a heater reaching the far corner of a room is also because the "carrier," air, transports heat via convection. Here's the interesting part—this term contains "velocity × velocity," making it nonlinear. That is, as the flow becomes faster, this term rapidly strengthens, making control difficult. This is the root cause of turbulence. A common misconception: "Convection and conduction are similar things" → They are completely different! Convection is carried by flow, conduction is transmitted by molecules. There is an order of magnitude difference in efficiency.
- Diffusion Term $\nabla \cdot (\Gamma \nabla \phi)$: Have you ever put milk in coffee and left it? Even without stirring, after a while, it naturally mixes. That's molecular diffusion. Now, next question—honey and water, which flows more easily? Obviously water, right? Honey has high viscosity ($\mu$), so it flows poorly. When viscosity is large, the diffusion term becomes strong, and the fluid moves in a "thick" manner. In low Reynolds number flows (slow, viscous), diffusion is dominant. Conversely, in high Re number flows, convection overwhelms, and diffusion plays a supporting role.
- Pressure Term $-\nabla p$: When you push the plunger of a syringe, liquid shoots out forcefully from the needle tip, right? Why? Because the piston side is high pressure, the needle tip is low pressure—this pressure difference becomes the force pushing the fluid. Dam discharge works on the same principle. On a weather map, where isobars are tightly packed? That's right, strong winds blow. "Flow arises where there is a pressure difference"—this is the physical meaning of the pressure term in the Navier-Stokes equations. A point of confusion here: "Pressure" in CFD is often gauge pressure, not absolute pressure. When switching to compressible analysis, if results suddenly become strange, it might be due to confusion between absolute/gauge pressure.
- Source Term $S_\phi$: Warmed air rises—why? Because it becomes lighter (lower density) than its surroundings, so it is pushed up by buoyancy. This buoyancy is added to the equation as a source term. Other examples: chemical reaction heat generated by a gas stove flame, Lorentz force applied to molten metal by an electromagnetic pump in a factory... These are all actions that "inject energy or force into the fluid from the outside," expressed by source terms. What happens if you forget a source term? In natural convection analysis, forgetting buoyancy means the fluid doesn't move at all—a physically impossible result where warm air doesn't rise in a winter room even with the heater on.
Assumptions and Applicability Limits
- Continuum Assumption: Valid for Knudsen number Kn < 0.01 (mean free path ≪ characteristic length)
- Newtonian Fluid Assumption: Linear relationship between shear stress and strain rate (non-Newtonian fluids require viscosity models)
- Incompressible Assumption (for Ma < 0.3): Treat density as constant. For Mach number 0.3 or higher, consider compressibility effects
- Boussinesq Approximation (Natural Convection): Consider density changes only in the buoyancy term, using constant density in other terms
- Non-applicable Cases: Rarefied gas (Kn > 0.1), supersonic/hypersonic flow (requires shock capturing), free surface flow (requires VOF/Level Set, etc.)
Dimensional Analysis and Unit Systems
| Variable | SI Unit | Notes / Conversion Memo |
|---|---|---|
| Velocity $u$ | m/s | When converting from volumetric flow rate for inlet conditions, pay attention to cross-sectional area units |
| Pressure $p$ | Pa | Distinguish between gauge and absolute pressure. Use absolute pressure for compressible analysis |
| Density $\rho$ | kg/m³ | Air: approx. 1.225 kg/m³ @20°C, Water: approx. 998 kg/m³ @20°C |
| Viscosity Coefficient $\mu$ | Pa·s | Note confusion with kinematic viscosity coefficient $\nu = \mu/\rho$ [m²/s] |
| Reynolds Number $Re$ | Dimensionless | $Re = \rho u L / \mu$. Indicator for laminar/turbulent transition |
| CFL Number | Dimensionless | $CFL = u \Delta t / \Delta x$. Directly related to time step stability |
Numerical Methods and Implementation
Details of Numerical Methods
Please tell me the numerical key points of the Eulerian granular model.
The biggest difficulty is handling when the solid volume fraction approaches the maximum packing fraction $\alpha_{s,max}$. In this region, solid phase pressure increases sharply, making it prone to numerical instability.
The Frictional Stress model is important; when $\alpha_s > \alpha_{s,min}$ (typically 0.5), frictional pressure and frictional viscosity are added using the Schaffer model or Johnson & Jackson model.
Drag Model Selection
How should I choose a drag model for gas-solid two-phase flow?
The Gidaspow model is the most common. It switches between the Ergun equation (dense packing region) and the Wen-Yu equation (dilute region) at a volume fraction of 0.8.
| Region | $\alpha_g$ | Model | Equation | ||
|---|---|---|---|---|---|
| Dense Packing | $< 0.8$ | Ergun | $\beta = 150 \frac{\alpha_s^2 \mu_g}{\alpha_g d_s^2} + 1.75 \frac{\alpha_s \rho_g \ | \mathbf{u}_g - \mathbf{u}_s\ | }{d_s}$ |
| Dilute | $\geq 0.8$ | Wen-Yu | $\beta = \frac{3}{4} C_D \frac{\alpha_s \alpha_g \rho_g \ | \mathbf{u}_g - \mathbf{u}_s\ | }{d_s} \alpha_g^{-2.65}$ |
Doesn't the discontinuity at the switch point cause problems?
Exactly, the discontinuity at the switch point can cause void fraction oscillations. The Huilin-Gidaspow model improves this by introducing a smooth transition function. The Syamlal-O'Brien model also uses a continuous equation across the entire region, offering higher stability.
Implementation in OpenFOAM
Which solver is used in OpenFOAM?
multiphaseEulerFoam supports Eulerian Granular. Each constitutive relation of KTGF can be selected in the kineticTheoryModel class. Main settings are done in constant/phaseProperties.
Settings in Fluent
In Ansys Fluent, enable Granular Phase within the Eulerian Multiphase Model. Important configuration items are as follows.
| Setting | Recommendation | Remarks |
|---|---|---|
| Granular Viscosity | Gidaspow | Standard |
| Granular Bulk Viscosity | Lun et al. | Bulk viscosity |
| Frictional Viscosity | Schaeffer | Dense packing region |
| Packing Limit | 0.63 | Random packing fraction |
| Restitution Coefficient | 0.9 | Typical value for glass beads |
The Wall of KTGF Convergence—The Battle with the "Negative Pressure" Phenomenon
The first stumbling block when implementing an Eulerian granular model is the problem where solid phase pressure becomes negative and diverges. This is a numerical instability that occurs when granular temperature rapidly converges to 0, making setting a floor value (minimum around 1e-10 m²/s²) a de facto essential countermeasure. ANSYS documentation presents the choice of "solving as a Partial Differential Equation or using an Algebraic approximation," but for high-density packing (α_s > 0.4), accuracy often cannot be achieved without the PDE method, requiring on-site judgment of the trade-off between computational cost and accuracy.
Upwind Differencing (Upwind)
1st order Upwind: Large numerical diffusion but stable. 2nd order Upwind: Improved accuracy but risk of oscillations. Essential for high Reynolds number flows.
Central Differencing (Central Differencing)
2nd order accuracy, but numerical oscillations occur for Pe number > 2. Suitable for low Reynolds number diffusion-dominated flows.
TVD Schemes (MUSCL, QUICK, etc.)
Suppress numerical oscillations while maintaining high accuracy using limiter functions. Effective for capturing shock waves or steep gradients.
Finite Volume Method vs Finite Element Method
FVM: Naturally satisfies conservation laws. Mainstream in CFD. FEM: Advantageous for complex shapes and multiphysics. Mesh-free methods like SPH are also developing.
CFL Condition (Courant Number)
Explicit method: CFL ≤ 1 is the stability condition. Implicit method: Stable even for CFL > 1, but affects accuracy and iteration count. LES: CFL ≈ 1 recommended. Physical meaning: Information should not travel more than one cell per time step.
Residual Monitoring
Convergence is judged when residuals for continuity, momentum, and energy equations drop by 3-4 orders of magnitude. The mass conservation residual is particularly important.
Relaxation Factors
Pressure: 0.2-0.3, Velocity: 0.5-0.7 are typical initial values. If diverging, lower the relaxation factor. After convergence, increase to accelerate.
Internal Iterations for Unsteady Calculations
Iterate within each time step until a steady solution converges. Internal iteration count: 5-20 times is a guideline. Residuals should...
Related Topics
なった
詳しく
報告