Incompressible Navier–Stokes, decomposed
RowSim solves the incompressible Navier–Stokes equations for a two-dimensional velocity field using Jos Stam’s stable-fluids decomposition, implemented as a sequence of GPU passes per frame. The governing system couples material acceleration, pressure gradient, viscosity, and body forces under the constraint that divergence remain zero (∇·u = 0). In practice the solver never attacks that system as a single monolithic kernel; it decomposes the update into stages that a real-time Metal pipeline can schedule, hazard, and branch by architecture.
Each frame applies, in order: semi-Lagrangian advection of velocity and density; injection of interaction and oar-trail forces; ripple injection, wave propagation, and ripple-to-velocity coupling; vorticity estimation and confinement; divergence assembly; pressure projection; and gradient subtraction. The projection step solves the Poisson equation ∇²p = ∇ · ũ with either a Red–Black Gauss–Seidel solver on Apple Silicon or a Jacobi solver on the fragment fallback. Without projection, injected energy becomes a smear. With it, the field stays continuous under room-scale projection.
That decomposition is the governing idea of the chapter: water is not a texture. It is a field that is moved, disturbed, measured for compression, and then corrected until the wake underfoot remains coherent from stroke to stroke. The rest of the fluid chapter makes that inventory of passes legible—advection’s stability bargain, vorticity’s restorative craft, pressure’s ownership model, ripples’ second memory, and the hybrids that keep floor projection readable.
Momentum
Incompressible Navier–Stokes
Momentum, divergence-free flow, and the pressure Poisson solve — the three relations the fluid chapter depends on.
Trace backward to move forward
Semi-Lagrangian advection asks where each texel came from. advect_compute traces each texel backward along the local velocity plus the constant carrier flow baseVelocity—default −2.0, right-to-left—then bilinearly samples the source field at the trace point. The bilinear sampler excludes taps that fall inside user-drawn obstacles through bilerpComputeAvoidUserObstacle, preventing density from leaking through solids that the obstacle mask has already declared off-limits to the wake.
Because the backward trace is unconditionally stable, the method permits large per-frame displacements without numerical blow-up—the defining property of Stam’s scheme—at the cost of implicit low-pass filtering. Numerical diffusion is the price of reliability at installation frame rates. Vorticity confinement and bloom stages partly compensate for that loss of sharpness perceptually, which is why the fluid chapter refuses to treat advection as a pure physics win or a pure aesthetic loss: it is the stability bargain that makes the rest of the pipeline possible.
Several non-physical stabilisers are layered on the advection result: a centre-distance dissipation term, an edge-bleed dampening when the carrier flow is near zero, and a flow-aligned hash-noise term gated by flowVisualizationEnabled that adds visible structure to an otherwise-uniform carrier. These are deliberate artistic–physical hybrids documented explicitly in the shader source. They do not claim to be viscosity; they claim to keep a right-to-left conveyor of water from looking empty under a projector.
Velocity field
How a 2D flow looks on a grid.
A synthetic curl-noise field for teaching — not the Metal advection or vorticity kernels. Click a mode to focus on flow, swirl, or speed.
01 · FlowEach arrow is a velocity sample. Direction is where the water moves; longer means faster.
Teaching velocity field (curl noise). Arrows show flow; tint shows swirl. Not the Metal vorticity kernel.
Inject energy, then restore curl
During a power stroke, anisotropic Gaussian splats inject velocity, density, and bloom energy at interaction and oar positions. Displacement-derived trail points are batched five per dispatch and splatted into the same fields, which is how the host’s trail geometry becomes a field disturbance rather than a polyline drawn on top of the water. Force magnitudes are capped per texel—maxForce approximately 7.2 to 7.8, with a global 15.0 cap in the force kernels—to prevent solver instability from overlapping splats when multiple producers fire into the same neighbourhood.
vorticity_compute estimates scalar curl from the staggered neighbourhood as roughly one-half the centred difference of the velocity components. vorticityConfinement_compute then computes a confinement force from the normalised gradient of absolute vorticity and adds it back to the velocity to counteract the rotational damping introduced by semi-Lagrangian advection. Without that restoration, the wake would advect stably but lose the swirling structure that makes a stroke read as water rather than as a blur.
Two refinements absent from the textbook formulation matter for hull legibility. Obstacle-edge curl shedding converts the tangential component of velocity relative to the carrier flow into additional rotational forcing where a texel neighbours the boat or cutout boundary, producing visible vortex shedding at the hull. Inside the boat hull, velocity is forced to a reduced carrier value using slipScale 0.62 rather than zero—a visual approximation of a no-penetration, partial-slip boundary that keeps the silhouette from looking like a frozen cutout in an otherwise living field.
Pressure makes the field incompressible
As noted in the governing decomposition, projection enforces ∇·u = 0 by solving ∇²p = ∇ · ũ and subtracting the pressure gradient. The divergence kernel assembles the right-hand side with centred first differences, substituting the local velocity component where a neighbour is an obstacle—a Neumann-like treatment that keeps solids from injecting false expansion into the Poisson problem. The solve uses the standard five-point Laplacian with α = −1.0 and β = 4.0.
On Apple Silicon, jacobi_rb_compute performs Red–Black Gauss–Seidel in place on the pressure slab. Each dispatch processes only one colour parity selected by a parityBuffer bound at buffer index one and compared against (gid.x + gid.y) modulo two. The renderer issues a red pass followed by a black pass per iteration, and sixteen such red-and-black pairs run each frame. In-place checkerboard update halves memory traffic and, because each colour’s neighbours were written in the preceding pass, converges roughly twice as fast per iteration as Jacobi.
On the fallback branch, jacobi runs twenty Jacobi sweeps with ordinary ping-pong buffers, reading one side and writing the other. The ownership diagram contrasts the fixed pressure slab under RBGS with the Jacobi x_in to x_out exchange. Naming both paths matters: the fluid looks continuous on either machine, but the pressure ownership model and iteration count are architecture-conditioned facts, not interchangeable implementation details.
∇²p = ∇ · ũ, then u = ũ − ∇p
Pressure solve
Checkerboard updates, faster convergence.
Brightness is pressure on a 14×9 grid. Red and black cells take turns so neighbours stay fixed during each sweep — about 2× faster than Jacobi for the same sweep count.
01 · RedUpdate every red cell from its four black neighbours. Half the grid advances in one parallel pass.
Checkerboard pressure solve on a 14×9 grid. Red then black sweeps; sixteen pairs on Apple Silicon.
A second field for wave memory
A separate ripple texture stores current and previous wave height in its red and green channels and evolves under an explicit damped wave equation: hₜ₊₁ = (2−d)hₜ − (1−d)hₜ₋₁ + c²∇²hₜ. A five-point Laplacian and a CFL clamp c ← min(c, 0.707) keep the scheme stable on a unit grid, since stability requires c² < 0.5. The wave equation is also advected along the carrier flow, so ripples drift with the fluid instead of freezing as a screen-space overlay that ignores the wake’s direction of travel.
Boundary handling zeros obstacle interiors as fixed boundaries, reflects top and bottom, and treats left and right viewport edges as open—substituting zero height—so waves leave the domain rather than reflecting forever. Those open lateral edges matter in a floor projection where leftover reflections would otherwise stack into visual noise. addRipple_compute injects Gaussian splats at up to five source positions with a small asymmetric secondary lobe for visual interest. The host de-duplicates pending ripples within a two-second window and adds a small random offset when a new ripple lands within twenty texels of an existing one, which prevents perfectly overlapping impulses from reading as a single stuck ring.
rippleToVelocity_compute couples the wave field back into the fluid: the height gradient drives a velocity perturbation of the form −∇h · s, and the temporal change hₜ − hₜ₋₁ drives a radial push, capped at 6.0 times rippleFluidStrength on the compute path and 3.0 times on the fragment path. That coupling is the mechanism by which roll and azimuth perturbations at the oars propagate into the visible wake, giving trunk orientation a fluid consequence rather than a UI badge.
Wave profile
A 1D slice of the ripple field.
Same damped wave update as WaveCompute.metal. Open ends, two pulse sites, defaults c=0.4 and d=0.02.
01 · HeightSolid curve is wave height h(x,t) right now — the surface the ripples ride on.
1D damped wave matching WaveCompute.metal. Ghost = previous frame; pulses ≈ oar entry (c=0.4, d=0.02).
Physics for continuity, craft for legibility
The stabilisers named earlier in this chapter stay on top of the core solve. Bloom energy accumulates at force-injection points, advects with velocity, decays with default bloomDecay 0.975 (flooring values below 0.001), and is blurred through a half-resolution thirteen-tap Gaussian before glow compositing. Particles, off by default, are a separate slab advected leftward at baseVelocity, seeded at the right edge with hash noise, and decayed so foam dots convey flow speed.
These are documented artistic–physical hybrids—not hydrodynamic fidelity certificates. Physics provides continuity; craft provides legibility: confinement, bloom, and the stabilisers that keep the carrier from reading as an empty conveyor. Palette mapping exists in the source and was off in every study session. The claim is an inspectable working field.