# GPU diffusion convergence: history and rationale

**Status (2026-09-26):** the implemented contract — operator and effective right-hand side, method, exits,
`DiffusionSolveSettings`, reporting convention, budget, kernel rules and tests — is the
[Diffusion README](/src/PolyCfd.Gpu/Diffusion/README). This page keeps why the solve audits the stored field
against an effective-RHS residual, why it has a bounded colored Gauss-Seidel fallback and one shared total budget, the
enclosure evidence behind those decisions, and why SOL-02a replaced the weighted-Jacobi primary phase described below
with Chebyshev semi-iteration (same operator, tolerance, budgets and fallback).

## Why an audited fallback and one shared budget

The relative residual is normalized by the **effective** right-hand side, including prescribed normal values,
tangential wall forcing and immersed terms; normalizing by the stored velocity alone rejected wall-driven starts from
rest (see [physics](/docs/PHYSICS_GUIDE#diffusion-and-walls)). Weighted Jacobi (ω = 0.8) alone stalled above the unchanged
`1e-6` target on the prescribed-flow enclosure at `nu*dt/h²=6.4`, and the former static 200 and moving 400 sweep budgets
did not cover every valid scenario (below). The 2026-09-24 decision kept the first 400 sweeps on the existing
weighted Jacobi path, continued unresolved solves with colored Gauss-Seidel on the same row evaluator and boundary
rules, audited the **stored float32 velocity** against the original effective right-hand side, and capped both phases
together at 800 iterations for static and moving steps. The fallback changes the iteration, not the backward-Euler
operator or its tolerance, and adds no retained buffers. Effective scenario settings record the budget; a resumed
campaign requires matching new evidence, not reuse of the old build's acceptance. Physical and performance
references and all scientific limits remain unchanged. See the [work board](/docs/WORKBOARD) for current validation
status and [performance validation](/docs/PERFORMANCE_VALIDATION).

## Chebyshev primary phase — 2026-09-26 (SOL-02a)

Weighted Jacobi dominated the moving-fan step: on the 128³ fan a 2026-09-25 experiment measured
276 sweeps and 287 ms of a 345 ms step. Chebyshev semi-iteration on the same undamped Jacobi row,
with the rigorous bound ρ = S/(1 + S), needed 56 sweeps (59 ms; step 124 ms), passed the quick
regression checks and changed fields by 1–3e-6. The user approved it (2026-09-25) as the one primary
method. SOL-02a kept the operator, tolerance, budgets, colored fallback and bytes, added deterministic
double audits with one host read each and an exit for inputs already within the tolerance, and adopted
the experiment's validated stall rule; bound-based exits are the proposed follow-up (SOL-05). On the
SOL-02a build the short fan's median diffusion count is 56 (Jacobi 266), the enclosure step-3 solve
that needed 468 iterations converges in 76, and circular Couette at d = 115 records the rule's early
exit as a known limit ([Diffusion README](/src/PolyCfd.Gpu/Diffusion/README#limits),
[complete validation](/docs/COMPLETE_VALIDATION#chebyshev-diffusion--2026-09-26)).

## Choice of fallback

The initial residual-correction prototype solved the first enclosure step at
`nu*dt/h²=6.4` in 560 sweeps, with independent residual `9.53e-7`. It still
stagnated at step 3 near `1.09e-6`, including a 2,400-sweep diagnostic. Longer
correction accumulation did not fix that plateau. The colored fallback is simpler
and completed all four 20-step enclosure startups (centered n64 dt=1/64,
centered n128 dt=1/256 and dt=1/128, offset n128 dt=1/128). At the formerly failing
third centered n128 dt=1/128 step, 468 sweeps give an independent stored-U residual
`9.8695035e-7`, matching the GPU report and production field exactly. The 92 CUDA
diffusion tests pass, including independent wall solutions, odd periodic grids and
honest nonconvergence reporting. Two complete before/after timing pairs pass all
scientific checks. Their combined verification-call time changes by +3.55%; the
first pair changes by +6.63% and fails its shear timing limit, while the reversed
pair changes by +0.32% and passes. Isolated warmed diffusion shows no sustained
penalty on the two tested fast-path grids. All attempts and the separate repeated
short cases are retained in the [measurement report] (unavailable local reference: `../../openfoam/compare/ENCLOSURE_N128_DIAGNOSTICS.md`).
That initial full CUDA run stopped at its historical timing gate after all 1,519 tests and nine scientific workloads passed. The ordinary-path cleanup removes runtime output checks/captured audit locals; 92 tests and all four startups pass again, and the cleaned solver passes a complete matched timing pair. Full validation completed in `output/validation/diffusion-full-verified-20260924`: all 1,525 selected tests, 25 native cases, channel refinement/replay and field regressions pass. Overall regression remains 30/31 because historical timing still fails after one fresh retry. The unchanged historical solver also fails its old gate under current conditions; monitoring caught intermittent Windows workers in both timing sequences. A quiet performance check is still needed. Startup success is not settled enclosure accuracy. See [the complete result](/docs/COMPLETE_VALIDATION#diffusion-fallback-validation--2026-09-24).

Near-unit periodic stress fields still demonstrate the limits of float storage.
An independent FFT solution has residual `7.21e-15`, but nearest float rounding
raises it to `1.339e-6`. This is not a lower bound over all float fields. Those
stress tests assert that convergence matches the actual stored-field residual;
they do not require every finite-budget solve to converge. No above-target field
is accepted. [Independent analysis] (unavailable local reference: `../../openfoam/enclosure/analyze_float_residual.py`)
retains its input/operator and limitations.

## Observed enclosure startup limits — 2026-09-24

The former static 200-sweep budget does not cover every valid scenario. Exact-input
CUDA probes on the prescribed-flow enclosure identify both iteration-budget and
float-update limits. At `nu*dt/h² = 3.2`, n64 dt=1/64 and n128 dt=1/256 require
272 and 280 sweeps to reach the unchanged target; their independent returned-field
residuals are 9.43e-7 and 8.84e-7. The n128 dt=1/128 case has coefficient 6.4:
its independent residual stalls at 1.27329e-6 through 2,400 sweeps because the
remaining corrections round away when added to the float32 velocities. Double row
accumulation is already used; these data identify a limitation of the former update,
not a proof that every float32 solution must miss the target.

The failure predates the latest format pull. Earlier n64 dt=1/256 enclosure
screening has coefficient 0.8 and converges in 76 sweeps. The full verify-all catalog
excludes this separate enclosure campaign. Pressure multigrid does not address
this diffusion limit. At diagnosis time no production tolerance/cap or protocol change had been made;
the expensive reference campaign remains paused pending validated solver and campaign evidence. See the [independent evidence and reusable probe] (unavailable local reference: `../../openfoam/compare/ENCLOSURE_N128_DIAGNOSTICS.md`).
These are first-step diagnostics, not completed physical comparisons.