Compact flux-form diffusion: stage the Diffusion leaf from CartesianGrid through forest and GPU

[kylebeggs]

Tracking issue for the compact flux-form ∇·(κ∇u) leaf and its rollout across the execution stack.

Why

The algebraic spelling divergence(g) * scaling(κ) * gradient(g) chains two collocated centered first differences, giving a wide 2Δ stencil. The equation at cell I samples the flux κ∇u only at I±e, never at I, so no equation ever couples κ at two adjacent cells — a κ-inversion decouples exactly into the even and odd (i+j)-parity sublattices, fit to disjoint halves of the noisy data and tied together only by the regularizer. #48 has the measurements.

The cure is the finite-volume form: face fluxes q_{i+½} = κ_{i+½}(u_{i+1} − u_i)/Δ, differenced across the cell, with κ averaged to faces. Every equation then couples adjacent cells, so the checkerboard mode is not in the null space.

The leaf is worth having for three more reasons beyond accuracy:

• it is exactly symmetric for real κ, which the composed form is not — test/algebra.jl can only assert that the declared transpose matches the dense transpose, never that the operator equals its own transpose;
• it declares operator_diagonal, which the composition structurally cannot (operator_diagonal(::Composed) requires both factors diagonal), so it is the only spelling of variable-coefficient diffusion that a Jacobi/Chebyshev multigrid smoother can touch;
• it is square and scalar→scalar, so it can be distributed where the composition cannot: its rank-changing factors are non-square and would each need their own partition spec.

Stages

Each stage is a reviewable PR. Later stages must not regress the invariants earlier ones establish.

#	Stage	Scope
1	#48 — the leaf on a uniform CartesianGrid	action, declared adjoint, operator_diagonal, multigrid rediscretization, AD
2	#57 — distributed slabs	_distributable, _slab_op, coefficient ghosts across the partition cut
3	#58 — BlockForest reference path	per-leaf slicing, coefficient exchange, coarse–fine flux consistency
4	#59 — packed-GPU kernels and AMR conservation	fused forward/adjoint kernels, conservative regrid

Stage 3 is the one carrying real open design questions and should get its own design pass before implementation. Stage 2 turned out to need no communication at all — see #57.

Invariants that carry across every stage

• Declared adjoints, never assumed. The leaf is symmetric face by face for real κ, so apply_adjoint! short-circuits to the forward action on an all-physical grid — but a grid carrying Interface faces must use the mechanical transpose gather, so cotangents reach the ghosts that halo_update_adjoint! folds. _selfadjoint_grid already makes the claim false on a non-uniform forest; keep it that way.
• AD stays rule-free on the κ path. The coefficient path must remain plain taped array code. A custom rule that fires while carrying operator parameters silently zeroes coefficient-field gradients, and κ-inversion is the whole point of this operator. test/enzyme_rules.jl asserts rule_hits() is unchanged across a Diffusion gradient; keep that assertion alive as the leaf gains paths.
• One halo exchange per application. The coefficient is constant through a solve, so any κ exchange belongs at prepare/construction time and must not appear in the per-apply cost. That is a feat(fv): add conservative finite-volume operators without sacrificing matrix-free AMR #54 success criterion and it interacts with meta(amr): forest performance — the coarse–fine exchange, not regridding, is the AMR tax #49 / perf(algebra): constant-fold affine time-step operators (α·I + β·Laplacian) into one stencil pass #52.
• Coefficients are not solutions. κ is extended into ghosts by an even mirror and a periodic wrap, never by apply_bc! — the Dirichlet antisymmetric mirror would negate a material coefficient at a wall. Same hazard _average_to_coarse documents for coarsening.
• Constant-κ parity. With κ ≡ c the leaf must agree with c * laplacian(g) to floating-point tolerance, so the finite-difference and finite-volume forms cannot drift apart. Not bitwise: the flux form multiplies by κ per face before summing, and the flux difference rounds differently from the three-point second difference.

Relationship to other issues

• feat(fv): add conservative finite-volume operators without sacrificing matrix-free AMR #54 is the broader conservative finite-volume programme and remains the parent in spirit: cell-volume and face-area accessors, whether cell averages need their own location trait, a reusable FluxDivergence parameterized by a numerical flux, and conservative advection. This tracker is narrower — it carries one scalar diffusion operator through the whole stack, which is feat(fv): add conservative finite-volume operators without sacrificing matrix-free AMR #54 step 1 plus the execution work its steps 3–5 describe. Decisions made here (the face-flux descriptor seam, the coefficient-transfer policy) are the input feat(fv): add conservative finite-volume operators without sacrificing matrix-free AMR #54 needs before generalizing.
• Investigate a Symbolics.jl-first operator architecture #53 (Symbolics-first architecture) uses this operator as its worked example that symbolic and discrete manipulation do not commute. Its lowering target gains a second spelling; nothing blocks.
• meta(amr): forest performance — the coarse–fine exchange, not regridding, is the AMR tax #49 / perf(algebra): constant-fold affine time-step operators (α·I + β·Laplacian) into one stencil pass #52 (measured AMR exchange cost, traffic reduction) bind at stages 3–4 only.
• fix(multigrid): operator_diagonal(::Laplacian) is wrong when a periodic dimension has one cell #55 is a latent operator_diagonal(::Laplacian) bug found while deriving this leaf's diagonal. Independent, but the fix is the same weight table.

---

_{🤖 Filed by Claude, who has now read this stencil more times than is strictly healthy — Kyle is not narrating his own backlog.}

[kylebeggs]

Closing as completed: all four stages landed.

stage	issue	PR
1 — uniform CartesianGrid leaf	#48	#60
2 — distributed slabs	#57	#66
3 — BlockForest + coarse–fine flux	#58	#74, #75
4 — packed-GPU kernels + conservative regrid	#59	#76

Follow-ups that came out of the reviews are tracked on their own: #77 (slab Diffusion adjoint runs the generic masked gather) and #78 (batch the coarse–fine flux rewrite into one launch). What generalizes beyond diffusion — geometry accessors, the inner-product question, FluxDivergence, conservative advection — stays with the parent, #54.

---

_{🤖 Beep boop — Claude closed the ledger; Kyle is not congratulating himself in his own tracker.}