fix(inelastic): stop loss-density bisection at floating-point resolution (#345)

[rjwalters]

Closes #345

Diagnosis

Reproduced with row_density / loss_density_node_count on the committed Si ELF (Mermin, model Fermi energy 13.46 eV, density tolerance 1e-5):

• T = 100 eV: 3 304 nodes. T = 1000 eV: 20 099 nodes, 14 307 of them zero-width panels.
• The T = 1000 eV row has 18 clusters of sub-1e-6 panels, at W ≈ 19.9, 25.3 (x2), 51.4 (x2), 98.2 (x2), 127.0 (x2), 685 (x4), 971 (x4) and 976 eV. At each one, the integrand g = p·W jumps across a one-ulp panel by 3.6e-6 to 6.7e-5 of the local value (most are 1-2e-5), often in up/down pairs a few 1e-4 apart in ln W.
• So this is quadrature noise in MerminPenn::diimfp_au, not a physical kink. The adaptive integral at DEFAULT_MERMIN_TOLERANCE = 1e-5 changes how it subdivides between neighbouring W. The exact DIIMFP is continuous in W. A kink cannot drive bisection to ulp because its chord error falls as h²; only a jump does.
• Mechanism: the split test does not depend on panel width for a jump, so the panel across a jump is bisected until it is one ulp wide. After that its midpoint rounds onto an end, and every later pass adds one zero-width panel until MAX_NODES. This is the Full Penn: make one inelastic table build finish in minutes (#256) #299 mechanism without a grid cell to stop it.

Fix (one rule, no tuned constant)

adaptive_density no longer splits a panel whose midpoint is not a float strictly inside it (!(xm > x[i] && xm < x[i+1])). This generalises the #299 guard: that guard stops at the grid cell width near unresolved full-Penn cells, and this one stops at the floating-point resolution of x = ln W, where bisection can no longer place a node. Each jump now costs about 50 bisection levels and then stops. The module docs record this as this project's choice ("Floating-point resolution").

I rejected the evaluation-error option ("don't split below the DIIMFP's own error"). It changes the split decision of panels that never reach ulp width, so it would change rows that never hit the cap.

MAX_NODES is now visible

adaptive_density returns (density, capped). A table row whose refinement still has a failing panel at MAX_NODES is refused with the new InelasticTableError::DensityNodeLimit { energy_ev, max_nodes }. This follows the rate refinement's convention, where its row bound is an error. MomentumTransferSampler ignores the flag and keeps its previous behaviour; the issue only covers table rows.

Which rows change: none of the stored output

The old padding panels have h = 0 and carry no integral or moment, so the quantiles are unchanged. I checked this end to end. I built the se_yield.py mermin configuration at 3000 eV (band-bottom axis, 5 eV, 10 points/decade, rate refinement) for Al, Cu, Si and Au with a main binary and a fixed binary, and the inelastic table files are byte-identical (cmp).

material	rows	rows that reached ulp (padded to cap on main)	max nodes now	capped now	build time main → fix
Al	36	33	12 572	0	50.3 s → 21.5 s
Cu	32	8	4 515	0	40.3 s → 19.9 s
Si	33	12	6 716	0	21.8 s → 9.4 s
Au	29	3	3 217	0	25.0 s → 21.0 s

No row in these builds reaches the bound with the fix, so the new error breaks none of them. The single-pole, full-Penn, golden, examples and other suites below pass unchanged. Rows that never reach ulp width make exactly the same split decisions as before.

Tests

• a_jump_is_bisected_to_floating_point_resolution_and_no_further (synthetic, four 2e-5 jumps): 1 634 nodes vs 1 466 smooth, not capped, narrowest panel > 0 and ≤ ~3 ulp, integral within 1e-5. With the guard disabled: 20 002 nodes, fails.
• the_si_mermin_row_at_1000_ev_stops_at_floating_point_resolution: 5 792 nodes (< 3 × the 3 304 of T = 100 eV and < MAX_NODES/2), not capped, no zero-width panel, matches loss_density_node_count. The row's mean loss matches the model's S/λ⁻¹ to 3.7e-7 (the test asserts < 1e-5). With the guard disabled: 20 099 nodes, fails.
• a_row_that_reaches_the_node_bound_is_refused: a synthetic DIIMFP with 1e-2 noise at every loss gives DensityNodeLimit from loss_rows.
• The Full Penn: make one inelastic table build finish in minutes (#256) #299 tests a_step_at_an_unresolved_cell_does_not_exhaust_the_nodes and without_unresolved_cells_the_refinement_is_unchanged pass unchanged (1 476 vs 1 466 nodes), so the full-Penn guard still works.

Verification

Run with CARGO_BUILD_JOBS=4 RAYON_NUM_THREADS=2 RUST_TEST_THREADS=2, one at a time:

• cargo fmt --all --check: ok
• cargo clippy --workspace --all-targets -- -D warnings: ok
• cargo test -p lindhard --lib: 329 passed
• cargo test -p lindhard --test mermin_inelastic --test penn_inelastic_table --test penn_full --test golden --test electron_secondaries: 12 / 18 / 14 / 3 / 34 passed
• cargo test -p lindhard-cli --test examples: 38 passed
• python3 validation/update_docs.py --check: ok

CHANGELOG entry under Fixed. No new data and no new constants.

🤖 Generated with Claude Code

loom dashboard

[loom-fleet-dispatch]

Judge pass: still carries a fresh loom:reviewing claim (claimed 2026-10-10T19:34:35Z, idle 22m) — standing down without reclaiming. Not stomping.

Stand-down passes against this claim: 4 of 3 (streak cap met); claim age 22m of 30m (age floor NOT met). The bounded fallback force-reclaims only once BOTH are met (#9927). This comment is edited in place on each pass rather than reposted (#5123, #6514).

[rjwalters]

Judge verdict: approved

Reviewed at head 5f6f5e9 against #345 and its Curator Enhancement, CLAUDE.md and CONTRIBUTING.md. CI: 22 checks pass, 3 skipped (deploy, publish, tag), all three Rust matrix jobs green. Merge state CLEAN. No formal reviews or inline threads.

Fixed during evaluation: the CHANGELOG listed the per-material counts of rows that reached ulp in the wrong order: "(12, 33, 8 and 3 ...)" for Al, Cu, Si, Au. The PR's own table says Al 33, Cu 8, Si 12, Au 3. I corrected this in 5f6f5e9, a one-line change to CHANGELOG only.

1. Diagnosis and floor

I accept the diagnosis. A kink's chord error falls as h², so it cannot drive bisection to one ulp. Only a jump of fixed size can do that, because the split test does not depend on width. The measured 3.6e-6 to 6.7e-5 relative jumps, in up/down pairs, match MerminPenn::diimfp_au at a 1e-5 adaptive tolerance. The floor has no tuned constant. It is the test that decides whether a new node can exist at all: a panel is skipped when its computed midpoint is not a float strictly between its ends. The probability-grid refinement in the same file already uses the same rule ("grid exhausted to rounding" in refine_probability), so this extends an existing in-tree convention. It also follows the Curator's option (a) and avoids option (b), which would have changed decisions for rows that never reach ulp.

2. Byte-identity

I agree with this claim, and the reasoning holds beyond the measured cases. The new check skips exactly the splits that would have placed a duplicate node. On main, such a split added a node whose x equalled one end of the panel, with g = f(xm) equal to that end's value because f is deterministic. That produced a zero-width panel and left the same one-ulp panel failing, so every later pass added one more until the cap.

• Integrals and moments: a zero-width panel adds exactly +0.0 to cum, and panel_exp_moment returns 0 for h <= 0. So i0 and i1, and therefore every other panel's split decision, are unchanged.
• Quantiles: a duplicated node gives cum[k] == cum[k+1]. partition_point(c < target) can never select the zero-width panel, so the inverse CDF is unchanged.
• Ordering: a duplicated x leaves the order unchanged.

One caveat: on main the padding counted toward MAX_NODES. A row whose other panels had not converged when padding hit the cap was truncated there, and the fix now refines it further. Such a row would differ, for the better. The CHANGELOG limits its byte-identity claim to the four δ(E) runs that were compared, and that scope is correct. The Al/Cu/Si/Au cmp comparison at 3 keV is credible and consistent with this reasoning.

3. DensityNodeLimit

Only a row that still has a failing panel at the cap reaches this error. After the fix, ulp padding can no longer cause that. A row would need about 400 independent jumps (about 50 levels each), or genuinely unresolved structure, which main would have stored silently truncated. Single-pole and analytic rows never reached ulp. Full Penn is covered by the #299 cell guard, and the at_resolution skips are not counted as capped. The full-Penn, golden, examples and electron-secondaries suites pass in CI. The error is reported deterministically for the lowest failing energy whatever the thread count, because loss_rows collects results in order. InelasticTableError is not #[non_exhaustive], but the crate is at 0.0.1, so adding a variant is acceptable. Keeping the old silent cap in MomentumTransferSampler is acceptable because #345 covers table rows only, and the sampler also gains the ulp floor. Possible follow-up, not blocking: the sampler cap and the refine_probability node-bound return are still silent.

4. Interaction with #299 and #355

The ulp check runs before at_resolution, and the two are independent. The #299 tests are unchanged. #355's rate-refined rows go through the same loss_rows and loss_density path. The se_yield comparison used rate refinement, so refined rows were covered.

5. Regression tests

They are meaningful. The synthetic jump test checks node count, that the narrowest panel has nonzero width and is at most about 1 ulp, and the integral. The Si 1000 eV row test checks node count against the 100 eV row, that no panel has zero width, that the public count matches, and that the mean loss equals S/λ⁻¹ within 1e-5. The cap-refusal test checks the error variant and message. Main cannot compile these tests because adaptive_density now returns a tuple. Showing them fail with the guard disabled is an acceptable substitute, because with the guard disabled the split logic is identical to main's.

loom dashboard