Ergun packed-bed pressure drop (Darcy -> Forchheimer)

Why this case matters

Nassu’s volumetric porous region applies a momentum sink with a linear (Darcy) coefficient porous_alpha and a quadratic (Forchheimer / canopy) coefficient porous_beta. The two earlier porous cases each exercise only part of this: case 02_porous_pipe_flow validates the linear Darcy term against the analytic Brinkman Bessel profile (porous_beta = 0), and case 01_flow_through_trees uses porous_beta on a canopy geometry with no clean analytic target. This case is the isolated pressure-drop sweep that activates both terms together and drives the flow across the linear -> quadratic regime transition, validating the combined sink against the Ergun packed-bed correlation Ergun[1].

It also carries a multiblock variant whose porous region straddles a refinement interface, guarding the issue #704 fix that restored the porous body force to the block-border population reconstruction.

The Ergun balance

A packed bed of spheres of diameter d_p and voidage eps (void fraction) develops a streamwise pressure gradient that the Ergun correlation splits into a viscous (Darcy) and an inertial (Forchheimer) contribution:

(1)\[\frac{dP}{dx} = 150\,\frac{(1-\varepsilon)^2}{\varepsilon^3}\,\frac{\mu\,u}{d_p^2} \;+\; 1.75\,\frac{(1-\varepsilon)}{\varepsilon^3}\,\frac{\rho\,u^2}{d_p}\]

with u the superficial velocity. In Nassu a volumetric region adds the body force \(F_\alpha = -\alpha\,u_\alpha - \beta\,|u|\,u_\alpha\) to the Guo source term on every region node. In a fully-developed, gradient-free porous-filled duct driven by a constant streamwise body force G, the viscous and advective terms vanish and the steady momentum balance is purely local:

(2)\[G = \alpha\,u + \beta\,u^2\]

This is (1) with \(dP/dx = \rho\,G\) (\(\rho = 1\)) once the lattice coefficients are identified (\(\mu = \rho\nu\), \(\rho = 1\)):

(3)\[\alpha = 150\,\frac{(1-\varepsilon)^2}{\varepsilon^3}\,\frac{\nu}{d_p^2}, \qquad \beta = 1.75\,\frac{(1-\varepsilon)}{\varepsilon^3}\,\frac{1}{d_p}\]

The two terms are equal at the transition velocity \(u^* = \alpha/\beta\), where the Forchheimer number \(Fo = \beta u/\alpha = 1\). Sweeping the drive G moves u from the Darcy-dominated regime (Fo << 1) through the transition to the Forchheimer-dominated regime (Fo >> 1).

Physical -> lattice mapping

Everything in the YAML is lattice units; the physical interpretation is recovered in the notebook from the parsed config.

Parameter

Lattice value

Note

cs^2

1/3

speed of sound squared

tau

0.56

relaxation time

nu

0.02

cs^2 (tau - 1/2)

eps

0.4

voidage (typical dense bed)

d_p

50

particle diameter (lattice nodes)

porous_alpha

6.75e-3

150 (1-eps)^2/eps^3 nu / d_p^2

porous_beta

0.328125

1.75 (1-eps)/eps^3 / d_p

u*

0.02057

transition velocity alpha/beta (Fo = 1)

ergunGate - the pressure-drop sweep

Setup

A 32^3 triply-periodic box, whole domain porous (D3Q27, RRBGK). A constant streamwise body force G drives the flow; with no gradients the steady state is the exact Ergun balance (2). The drive is swept via !unroll over F.x, with porous_alpha / porous_beta baked once and shared across the sweep (v1 bakes one coefficient pair per kernel).

Velocity / Reynolds sweep

Each G is \(\alpha u + \beta u^2\) solved for the target u; sim_id follows the table order.

sim_id

u_target

Fo

Ma

G = F.x

regime

0

0.004

0.194

0.0069

3.225e-05

Darcy-dominated

1

0.010

0.486

0.0173

1.00312e-04

linear, near-onset

2

0.020

0.972

0.0346

2.6625e-04

transition (Fo ~ 1)

3

0.035

1.701

0.0606

6.38203e-04

Forchheimer-leaning

4

0.050

2.431

0.0866

1.15781e-03

Forchheimer-dominated

All Ma < 0.1 and tau in (0.5, 2).

Validation metrics

For each sim_id the steady superficial velocity u is the domain-mean streamwise velocity at the last snapshot. The pair (u, dP/dx = G) must land on the analytic Ergun curve (1) to within ~3% across the whole sweep, with the curvature of dP/dx(u) (the quadratic Forchheimer term) resolved at the high-Fo end. A point probe records the exponential approach to the terminal velocity.

ergunRefinedSlab - multiblock interface guard

Setup

The same Ergun box driven at the transition point (u_target = 0.02, Fo ~ 1, so both terms are active), but with an inner slab in x (x in [8, 24]) refined to level 1. The porous region straddles two C2F / F2C refinement interfaces.

Validation metrics

This is the regression guard for issue #704 (the porous force was omitted from the block-border population reconstruction, producing a per-block sawtooth at interfaces). Because the Ergun balance is purely local and level-independent, the converged streamwise velocity must be spatially uniform - the same value on the coarse and fine blocks, with no sawtooth at x = 8 or x = 24, and matching the sim_id 2 value of the gate sweep within ~1%.

Note

The notebook comparison and the committed figures require a GPU run; results are pending. The configs parse and are physically well-posed.