AeroSim Docs
Porous-Media Momentum Sink¶
Most of the flows Nassu targets are bounded by solid walls and open boundaries, but two recurring situations call for a resistance that is spread through a volume rather than imposed on a surface. The first is a genuine porous medium: a vegetation canopy, a perforated screen, a packed bed, or any permeable obstruction that the flow passes through while losing momentum to the solid matrix. The second is purely numerical: a sponge layer placed just ahead of an outlet that absorbs outgoing pressure waves and velocity overshoots before they can reflect back into the domain and contaminate the solution.
Both needs are met by the same mechanism, a volumetric momentum sink applied to every fluid node inside a tagged region, combining a linear (Darcy) term and a quadratic (Forchheimer / canopy) term. This page sets out the physics that motivates the sink, the closed-form solutions used to validate it, and how it enters the lattice Boltzmann update. The key point is that it requires no new collision operator and no new stored quantity: it is a body force, and it rides the same Guo forcing route the solver already uses for the global driving force (see the discrete-velocity force).
Key idea: one mechanism, several uses
A region tagged as porous feels a body force that opposes its own velocity, with a part linear in the velocity, \(-\alpha\, u_\alpha\), and a part quadratic in the speed, \(-\beta\,|u|\,u_\alpha\). With a large \(\alpha\) the region behaves as a strongly resistive porous block or an outlet sponge; with a small \(\alpha\) it behaves as a weakly permeable medium that the flow threads through; the quadratic \(\beta\) term adds the speed-squared resistance of a vegetation canopy or a high-speed packed bed. The physics is identical across uses; only the magnitudes of \(\alpha\) and \(\beta\) distinguish a canopy from a wave-damping layer.
Darcy’s law and the linear momentum sink¶
Flow through a porous medium at low pore-scale Reynolds number obeys Darcy’s law: the volume-averaged velocity is proportional to the pressure gradient, with the constant of proportionality set by the permeability \(K\) of the medium [1]. Equivalently, the medium exerts a drag force on the fluid that is linear in the local velocity,
where \(\nu\) is the kinematic viscosity and \(\alpha = \nu / K\) is the linear Darcy resistance coefficient (lattice units, \(\alpha \ge 0\)). This is exactly the term Nassu adds. The resistance is a pure sink: it always opposes the motion, it removes momentum at a rate proportional to how fast the fluid is moving, and it carries no preferred direction of its own.
The strength of the medium is most naturally expressed through the dimensionless Darcy number, the permeability scaled by a characteristic length \(L\),
A large \(\mathrm{Da}\) means a weakly resistive, highly permeable medium that barely perturbs the flow; a small \(\mathrm{Da}\) means a strong resistance that nearly arrests it. The same parameter controls how an outlet sponge behaves: a strongly damped layer is just a small-\(\mathrm{Da}\) porous block.
Why a sink and not a fixed velocity
A momentum sink lets the flow find its own balance rather than imposing one. In an outlet sponge this matters: the layer must absorb whatever disturbance arrives without dictating the velocity, so that the interior solution is unconstrained. The linear-in-\(u\) form guarantees that a stronger disturbance is damped harder while a quiescent region is left untouched, which is precisely the behaviour a non-reflecting layer needs.
The Forchheimer (canopy) quadratic drag¶
Darcy’s law holds only while the pore-scale Reynolds number is small. As the through-flow speeds up, inertial losses in the pore network grow faster than linearly, and the resistance acquires a term proportional to the square of the velocity. The Forchheimer extension captures this by adding a quadratic drag to the linear Darcy sink,
where \(|u| = \sqrt{u_\beta u_\beta}\) is the local speed and \(\beta \ge 0\) is the quadratic (Forchheimer) drag coefficient (lattice units). The factor \(|u|\) makes the term grow with the square of the speed while keeping it aligned with, and opposed to, the velocity, so like the Darcy term it is a pure sink that carries no preferred direction of its own.
The same quadratic form is the natural model for canopy drag. A vegetation canopy resists the wind with a force per unit volume
set by a drag coefficient \(C_{\mathrm{D}}\) and a leaf (or frontal) area density \(a\) (area per unit volume), which maps onto Equation (3) with \(\beta = C_{\mathrm{D}}\, a\). This is why a single quadratic coefficient serves both a high-Reynolds packed bed and a forest canopy: the speed-squared dependence is the same, and only the magnitude of \(\beta\) distinguishes them.
Citation needed
The Forchheimer and canopy-drag forms above are stated generically; a specific reference for the quadratic packed-bed correction and the leaf-area-density canopy model should be added to the bibliography when one is settled on.
The Brinkman extension and the wall¶
Pure Darcy drag, Equation (1), has no viscous term, so on its own it cannot honour a no-slip wall: the velocity would jump to its terminal value right up to the boundary. The Brinkman extension repairs this by retaining the viscous diffusion term alongside the Darcy resistance [2], so that the steady, fully developed momentum balance under a constant driving force \(G\) reads
between the driving force, viscous diffusion, and the linear sink. In Nassu the Darcy term is the body force added in the region and the viscous term is already present in the LBM collision, so the Brinkman balance is recovered automatically wherever a porous region abuts a no-slip wall. The result is a thin viscous adjustment layer, the Brinkman boundary layer, whose thickness scales with \(\sqrt{\nu/\alpha}\) and across which the velocity drops from its porous-core value to zero at the wall.
This balance has two instructive limits. Far from any wall the viscous term is negligible, the driving force exactly cancels the Darcy resistance, and the velocity settles to the terminal value
In the opposite limit of vanishing resistance (\(\alpha \to 0\), \(\mathrm{Da} \to \infty\)) the sink disappears and the balance reduces to ordinary viscous Poiseuille flow.
How it enters the LBM¶
The porous sink is added to the lattice Boltzmann update as a body force through the Guo forcing scheme Guo et al.[3], the same route the solver uses for its global driving force \(F\) and for the buoyancy force. On a node flagged as porous the per-direction Cartesian force
is folded into the total body force, which the collision then expands onto the discrete velocities with the standard Guo discrete-velocity force (see Equation force on the LBE page). Because both terms are functions of the already-loaded velocity \(u_\alpha\), they need no extra macroscopic of their own: the bulk kernel reads the velocity it has already computed, forms the speed \(|u|\), multiplies by the baked-in \(-\alpha\) and \(-\beta\), and adds the result to the source term. The Guo half-step correction to the velocity (described on the macroscopics page) applies unchanged.
No new collision path
The porous region reuses the entire fluid collision unchanged. There is no porosity-modified equilibrium and no second collision operator; the medium is represented purely as an additional Guo body force gated on a per-node flag. The resistance coefficients \(\alpha\) and \(\beta\) are compile-time constants baked into the generated kernel, so a porous region is the same specialised kernel as the bulk fluid with one extra term, carrying no per-node storage and no runtime branch. Both coefficients are level-0 lattice values, rescaled by \(1/2^{\mathrm{lvl}}\) per refinement level so the physical resistance is independent of the local grid spacing.
Warning
The linear \(\alpha\) term models Darcy drag and is faithful in the Darcy and Darcy-Brinkman regimes; the quadratic \(\beta\) term adds the Forchheimer / canopy correction that becomes important at higher pore-scale Reynolds number. Keep \(\alpha \ge 0\) and \(\beta \ge 0\): a negative coefficient would inject momentum rather than absorb it.
Analytic validation anchors¶
Two exact solutions anchor the implementation, both reproduced by the Darcy-Brinkman porous pipe validation case.
Uniform Darcy balance. A periodic box filled entirely with porous medium and driven by a constant force \(G\) has no walls and no viscous term in the balance, so the steady state is the exact Darcy force balance \(G = \alpha\, u\). The velocity must approach the terminal value \(u = G/\alpha\) of Equation (6) everywhere, which directly checks that the sink is applied with the correct sign and magnitude.
Brinkman pipe profile. A circular pipe of radius \(R\) filled with porous medium and driven by a constant axial force \(G\), with no-slip at the wall, has the closed-form axial velocity profile that solves the Brinkman balance (5) in cylindrical coordinates,
where \(I_0\) is the modified Bessel function of the first kind, order zero. The shape is set by the Darcy number \(\mathrm{Da} = \nu/(\alpha R^2)\): at large \(\mathrm{Da}\) the Bessel term is small and the profile approaches the parabolic Poiseuille shape, while at small \(\mathrm{Da}\) the core flattens to the plug value \(G/\alpha\) and the velocity drops to zero only across a thin Brinkman layer at the wall. Sweeping \(\alpha\) across this range and comparing against Equation (8) confirms that the balance between driving force, Darcy resistance, and viscous wall drag is recovered.
The outlet-damping use case¶
The same volumetric sink doubles as a numerical sponge. Placing a porous slab over the last portion of the streamwise extent, just ahead of a zero-gradient outlet, damps the pressure waves and velocity overshoots that an impulsively started or strongly disturbed flow would otherwise drive into the open boundary. Without the slab those disturbances reflect off the outlet and travel back upstream; with it, the linear sink absorbs them in proportion to their amplitude, so the interior solution is left undisturbed.
Use case: a non-reflecting outlet layer
For an outlet sponge, choose \(\alpha\) large enough to damp the disturbance within the slab but graded gently enough that the entry to the porous region does not itself reflect. The slab is configured as a volumetric region with a pos predicate selecting the near-outlet band; the resistance is the same \(\alpha\) used for a physical medium, only larger. Because the term is a body force and nothing more, the sponge composes freely with the existing outlet boundary condition.