Crate pantometry_fluid 

Incompressible flow in three dimensions, by projection on a staggered grid.

  ∂u/∂t + ∇·(uu) = −∇p/ρ + ν∇²u + g
  ∇·u = 0

Why this is the hardest of this workspace’s domains to trust

Every other physics here has closed forms lying around. Fluids has few, its schemes trade stability against numerical diffusion, and “it looks like a fluid” is the easiest wrong answer in computational physics to accept: a scheme with the wrong viscosity still makes plausible vortices, and a scheme that quietly loses momentum still makes a pretty picture.

So this crate is built around the three exact solutions that exist, and each one is chosen to be blind to a different mistake:

  Poiseuille    u(y) = (g/2ν)·y(h−y)        exact — a quadratic, and a second difference of one
  Couette       u(y) = U·y/h                 exact — linear, and blind to advection entirely
  Taylor–Green  e^{−2νk²t}                   the full nonlinear equations, decay rate and all

The first two are unidirectional and steady, so the advection term is identically zero in both — they cannot check it at all, and saying so is the point. Taylor–Green can: it is an exact solution of the complete equations, in which the nonlinear term is balanced by the pressure gradient rather than absent. Beside it sit two statements that hold at machine precision and catch what a decay rate is too coarse to see: a uniform flow must remain exactly uniform, and total momentum in a periodic box must not move at all.

The staggered grid, and its one guarantee

Velocities on cell faces, pressure at cell centres — the same arrangement Yee uses for electromagnetism, and for the same reason: the divergence of a face velocity lands naturally at a cell centre, and the gradient of a centred pressure lands naturally on a face. No interpolation, and no checkerboard pressure mode.

After the projection, ∇·u is the residual of the pressure solve and nothing else — see Channel::divergence. That is weaker than electromagnetism’s identity, which holds exactly: here it holds to whatever the conjugate-gradient solve was asked for.

Two limits, and one of them is on the grid rather than on the step

  dt ≤ dx²/(6ν)      viscous, the same Fourier limit conduction has
  dt ≤ dx/|u|max     advective, the Courant limit
  |u|dx/ν ≤ 2        the cell Reynolds number — a property of the *mesh*

The third is the one that surprises people. Central differences on the advection term go unstable when a cell is too coarse for the viscosity to smooth what advection sharpens, and no amount of shortening the step fixes it: the mesh is wrong. Channel::cell_reynolds reports it and Channel::step refuses above two, rather than producing the wiggles that a reader would take for turbulence.

What is deliberately not here

No turbulence model, no compressibility, no free surface, no immersed geometry, no adaptive mesh. A box with periodic sides and optional walls, which is exactly what the three exact solutions live in.

Structs

Channel

A rectangular box of incompressible fluid.

Fluid

A Newtonian fluid.

Enums

Walls

What the y faces of the box are.

Constants

CELL_REYNOLDS_LIMIT

The largest cell Reynolds number central differences stay stable at.

Functions

poiseuille_mean_speed

The mean speed of plane Poiseuille flow driven by a body force.

taylor_green_rate

The rate at which a Taylor–Green vortex’s velocity decays, 2νk².