pushpull: export pull/push/count/grad_backward + a finite-difference oracle

[balbasty]

Part of the stack for fastfields-torch#16 (backpropagate through the pushpull sampling grid).

Stacked PR. Base is claude/cuda-pushpull-dynamic-spline (#48), not main — this uses its FF_BOUND_*/FF_SPLINE_* selector macros and BoundVec/SplineVec carriers. The impl pin points at fastfields-cpu-impl#39, also unmerged; both need re-bumping once the parents land.

What

The backward kernels have existed in the impl layer since the port but were never reachable — nothing above cpu-impl could call them, so fastfields.torch had to refuse to differentiate through grid.

This exports all four through the public ABI, mirroring how pull/push/count/grad are already exposed: same DISPATCH_PP matrix, same dtype/ndim/bound/order handling, same 32-bit index narrowing.

Contract (documented in pushpull.h):

out	d/d(inp). Accumulated into where the op scatters (pull_backward, grad_backward), overwritten where it gathers (push_backward). Callers should always pre-zero it.
gout	d/d(grid). Always overwritten — one grid point maps to exactly one output element.
ginp	the incoming gradient, shaped like the forward op's output.
out-of-FOV	contributes nothing, gets a zero gout.

abs is exposed only on grad_backward, where it must match the forward grad call. The pull/push/count adjoints pin it to false: abs swaps the signed spline derivative for its absolute value, a majorisation trick that only makes sense for the grad operator itself, and a true adjoint always needs the signed one. That also halves what would otherwise be a doubled instantiation matrix.

Build layout

The adjoints get their own translation unit, pushpull_backward.cpp, because they instantiate the whole ndim × order × bound × dtype matrix a second time — 7 min on its own at the library's fully-static policy. Worth making parallelisable rather than serialising it behind pushpull.cpp; same reasoning as the CUDA reg_field/reg_field_rls split.

The dispatch macros the two units must agree on are factored into a new private pushpull_dispatch.h rather than copy-pasted, so they cannot drift on which (order, bound) pairs are static vs. Dynamic. That is the one non-additive change in this PR: pushpull.cpp loses its macro block and gains an #include.

The oracle test

tests/test_pushpull_backward.cpp differences the scalar loss

L(inp, grid) = <forward_op(inp, grid), ginp>

against every element of both inp and grid, and compares to what the backward op returns. That is a full Jacobian check against the independently-tested forward ops, not a smoke test.

6381 checks over ndim 1/2/3, orders 0–3, bounds DCT2/DST2, with and without a batch dim — plus direct assertions for the out-of-FOV gate (a step function, which no finite difference can see), the accumulate-vs-overwrite semantics, and a transposed (non-contiguous) field.

Two deliberate choices worth flagging for review:

• Field and grid get different spatial shapes. Equal shapes make several stride arrays coincide, which is precisely how the push_backward stride bug survived in jitfields.
• Grid coordinates are kept off the spline knots. A spline of order k is only C^(k−1), so the higher derivatives genuinely jump at a knot and a central difference straddling one would be meaningless — a property of the interpolant, not of the implementation.

Bugs this found

Four, all inherited from jitfields, all fixed in the commits this pins:

1. push_backward indexed ginp with stride_inp — 801/6381 failures (cpu-impl#39)
2. 1D linear grad_backward returned a zero gradient wrt inp (kernels#45)
3. 2D linear grad_backward returned a zero gradient wrt grid (kernels#45)
4. 3D linear grad_backward — same (kernels#45)

Verification

make test: 12/12 binaries, 34306 checks, 0 failures.

Workstream: claude-jitfields-to-fastfields

🤖 Generated with Claude Code

---

Generated by Claude Code

[balbasty]

good to go