diff --git a/.claude/skills/nfft-perf-eng/REFERENCE.md b/.claude/skills/nfft-perf-eng/REFERENCE.md new file mode 100644 index 00000000..dc49ef72 --- /dev/null +++ b/.claude/skills/nfft-perf-eng/REFERENCE.md @@ -0,0 +1,104 @@ +# Performance optimization loop (test-pinned, benchmark-measured) + +Overview and map. Each phase has its own detail doc under [`details/`](details/) — +read the phase doc when you reach that phase ([SKILL.md](SKILL.md) is the checklist). + +How a coding agent optimizes a scoped code region (a function, a hot loop) in this +repository **without human intervention**, keeping two independent feedback signals +in hand the whole time: + +- **Correctness** — the CUnit accuracy suites (see + [`test-methodology.md`](../../../docs/agents/test-methodology.md)). +- **Performance** — the CodSpeed / google_benchmark binaries built by the **CMake** + tree, measured locally in **walltime** mode (see + [Measurement modes](details/measurement-modes.md)). + +The core idea: before changing anything for speed, *prove* which tests catch a +regression in the target and which benchmarks measure its speed. Do that by +deliberately breaking the target and watching what reacts. Then — still before any speed +edit — *analyze* the target's rounding error under the standard model, so the optimization +is held to the right accuracy (the derived bound, not merely a passing test) and an +*avoidable* error source is rooted out rather than preserved. Then you optimize with a +known, trustworthy net (tests), a known yardstick (benchmarks), and a known error model. + +The shape of the whole task: + +``` +Step 0 open the deliverables dir ── gitignored .perfeng/ + tracker README + squash base +Phase A full baseline ── build walltime trees; capture EVERY test + benchmark result (the exit reference) +Phase B pin correctness net ── which tests guard the target? [HARD GATE] +Phase C pin performance metric ── which benchmark measures it? [HARD GATE] +Phase D rounding-error analysis ── standard-model bound; set the accuracy objective +Phase E inner loop ── optimize toward that objective against the scoped net + metric (commit each unit) +Phase F exit gate ── re-run the FULL Phase-A baseline; no failure, no regression +Phase G conclude ── squash to one commit; offer push + PR (label perf-eng); package .perfeng-pr- zip + attach +``` + +**Front-to-back tracking.** The whole run lives in one fixed, **gitignored** directory, +`.perfeng/` — a tracker `README.md` plus one deliverable per phase. Deliverables are **never +committed**; the run's record ships as its squashed commit + PR + the zip attached to that PR +(Phase G). Each phase produces concrete, canonical-format deliverables, and *deliverables = exit +gate*: a phase isn't done until its deliverables are written and the tracker row is flipped. The +layout, tracker template, and snapshot formats are in +[`details/deliverables.md`](details/deliverables.md) — read it before Phase A. + +Phases B and C carry **hard gates**: if breaking the target fails *no* test, or if +*no* concrete benchmark metric can be obtained for it, the task **cannot proceed** — +that is a coverage gap to fix or escalate, not something to work around. The narrow +net/metric from B and C drive the fast inner loop (Phase E); the full Phase-A +baseline is the slow, authoritative check at the end (Phase F) that nothing outside +that narrow scope was broken or slowed down. Between the gates and the inner loop, +Phase D analyzes the target's rounding error and sets the *accuracy* objective the inner +loop optimizes toward (root out avoidable error first, or hold the derived bound). + +Even that authoritative check is *bounded* — it covers specific sizes and grids. So every +run also carries a **risk assessment** ([`details/risk-assessment.md`](details/risk-assessment.md)): +a survey of the accuracy side effects the tests in hand can't see (e.g. a precision loss +that only shows at larger transforms), proven or not, surfaced in the human summary. When a +risk is material and cheaply testable, the net is **extended** +([`details/extending-tests.md`](details/extending-tests.md)). + +The worked example throughout is `X(trafo_direct)()` — the direct, O(N·M) NDFT +(`kernel/nfft/nfft.c:145`): a real transformation, benchmarked forward and adjoint in +1d/2d/3d. Picking a *wrong* target is possible (code no test pins, or that no benchmark +measures), but you don't have to detect that up front — the Phase B and C hard gates +surface it for you: no failing test, or no benchmark moves, means **stop**. See +[caveats](details/caveats.md) for the common ways a target turns out unmeasurable. + +## Map — read these in order + +| Step | Detail doc | +|------|-----------| +| **Step 0 / deliverables** — task directory, tracker, canonical formats | [`details/deliverables.md`](details/deliverables.md) | +| **Phase A** — build tree + full baseline | [`details/phase-a-baseline.md`](details/phase-a-baseline.md) | +| **Phase B** — pin the correctness net *[HARD GATE]* | [`details/phase-b-correctness-net.md`](details/phase-b-correctness-net.md) | +| **Phase C** — pin the performance metric *[HARD GATE]* | [`details/phase-c-performance-metric.md`](details/phase-c-performance-metric.md) | +| **Phase D** — rounding-error analysis | [`details/phase-d-error-analysis.md`](details/phase-d-error-analysis.md) | +| **Phase E** — inner loop | [`details/phase-e-inner-loop.md`](details/phase-e-inner-loop.md) | +| **Phase F** — exit gate | [`details/phase-f-exit-gate.md`](details/phase-f-exit-gate.md) | +| **Phase G** — conclude (squash · push · PR · attach) | [`details/phase-g-conclude.md`](details/phase-g-conclude.md) | + +Cross-cutting references, consult as needed: + +| Topic | Doc | +|-------|-----| +| task directory layout, tracker, canonical deliverable formats | [`details/deliverables.md`](details/deliverables.md) | +| **float · double · long double** matrix (three trees; A/D/E run all) | [`details/precision-matrix.md`](details/precision-matrix.md) | +| fill-in skeletons to copy per phase (tracker + one per deliverable) | [`templates/`](templates/) | +| **walltime** is the metric (local + CI); **the noise rule**; the optional callgrind tie-breaker | [`details/measurement-modes.md`](details/measurement-modes.md) | +| pitfalls (crash faults, OpenMP-only changes, narrow benchmark coverage, benign errors) | [`details/caveats.md`](details/caveats.md) | +| **risk assessment** — what a green run can still hide; the risk table | [`details/risk-assessment.md`](details/risk-assessment.md) | +| **extending the net** — online + refgen test cases to prove/disprove a risk | [`details/extending-tests.md`](details/extending-tests.md) | +| what is agent-operable (the whole walltime loop, no account) | [`details/tooling-status.md`](details/tooling-status.md) | + +## Operational reference + +| Step | Command | +|------|---------| +| Configure (local loop) | `cmake -S . -B build-cmake -DNFFT_BENCHMARK_MODE=walltime -DNFFT_ENABLE_OPENMP=ON` | +| Build (lib + tests + benchmarks) | `cmake --build build-cmake -j` | +| Run tests | `ctest --test-dir build-cmake` | +| Granular test run | `build-cmake/tests/checkall` / `…/checkall_threads` (exit code + stdout `-> OK/FAIL`) | +| Failing cases only | `build-cmake/tests/checkall 2>&1 \| grep -E '\-> (FAIL\|ERROR)'` | +| Measure (walltime — the metric, local + CI) | `CODSPEED_PROFILE_FOLDER=/tmp/b build-cmake/benchmarks/bench_nfft_direct --benchmark_filter='…'` → `/tmp/b/results/*.json` (`median_ns`) | +| Tie-breaker only (layout vs real work, no account) | build `-DNFFT_BENCHMARK_MODE=simulation`; `valgrind --tool=callgrind … --benchmark_filter=''` → process-total `I refs` delta ([measurement-modes](details/measurement-modes.md)) | diff --git a/.claude/skills/nfft-perf-eng/SKILL.md b/.claude/skills/nfft-perf-eng/SKILL.md new file mode 100644 index 00000000..ec877065 --- /dev/null +++ b/.claude/skills/nfft-perf-eng/SKILL.md @@ -0,0 +1,183 @@ +--- +name: nfft-perf-eng +description: Safely optimize a scoped C region in the NFFT3 library for performance and/or accuracy without human intervention. Tests (CUnit) and benchmarks (CodSpeed walltime) ensure correctness and performance as independent signals. Deterministic steps (build, capture, compare, trend fit) are driven by helper scripts. Use when asked to optimize, speed up, or reduce the cost of a specific function/loop/kernel in this repo. +--- + +# Optimize NFFT performance (test-pinned, benchmark-measured) + +Optimize a scoped region with two independent signals in hand the whole time: +**correctness** (CUnit accuracy suites) and **performance** (CodSpeed/google_benchmark +binaries). Before changing anything for speed or accuracy, *prove* which tests catch +a regression and which benchmark measures the speed — by deliberately breaking the target and +watching what reacts. Then optimize with a trustworthy net and a known yardstick. + +**Full methodology:** [REFERENCE.md](REFERENCE.md) is the overview + map; each phase +below links to its own detail doc under [`details/`](details/). This SKILL is the +enforced checklist — open a phase's doc when you reach it, plus the cross-cutting +[deliverables](details/deliverables.md), [precision-matrix](details/precision-matrix.md), +[measurement-modes](details/measurement-modes.md), [caveats](details/caveats.md), +[risk-assessment](details/risk-assessment.md), [extending-tests](details/extending-tests.md) +and [tooling-status](details/tooling-status.md) docs as needed. + +**Every phase produces deliverables.** The loop is tracked front-to-back in one fixed, +**gitignored** directory, [`.perfeng/`](details/deliverables.md#where-it-lives) — **deliverables +are never committed.** The run's permanent record ships as its squashed commit + the PR + the +zip attached to that PR (Phase G), not as repo files. +*Deliverables = exit gate:* a phase is not done until its deliverables exist in `.perfeng/` +and the tracker row is updated. The format, layout, and canonical snapshot shapes are in +[deliverables.md](details/deliverables.md) — read it before Phase A. Don't hand-write +deliverables: each phase has a fill-in skeleton under [`templates/`](templates/); +copy the named one and fill its `<…>` placeholders. + +**Commit source changes as you go.** `.perfeng/` is gitignored, so commits capture only +*source* changes — commit each self-contained unit of work (a kept optimization step, a +permanent test addition) with a clear message as you reach it. Phase G squashes these into one +at the end. + +## The loop (do these in order) + +Create one TodoWrite item per phase. **Phases B and C are HARD GATES — if they fail, +stop and report, do not work around them** (record the gate failure as a `blocked` +deliverable — see [deliverables.md](details/deliverables.md)). + +- [ ] **Step 0 — open the deliverables directory.** Run + [`scripts/perf-init.sh `](scripts/perf-init.sh): it creates the fixed, gitignored + `.perfeng/`, copies every template in (tracker → `README.md`, the phase docs, + `error-analysis.html`), stamps the baseline commit, and records the **squash base** (HEAD now) + into `.perfeng/BASE` for Phase G. Set the **Target** line. → *Deliverable:* the `.perfeng/` + directory + initialized tracker. + +- [ ] **[Phase A — baseline](details/phase-a-baseline.md).** Build all three precision trees in + **walltime** mode and capture the FULL test + benchmark state (the exit reference for Phase F) — **in + all three precisions** (the sources are precision-agnostic; see [precision-matrix](details/precision-matrix.md)). + Both steps are deterministic: + ```bash + SCR=.claude/skills/nfft-perf-eng/scripts + $SCR/perf-build.sh walltime # configure+build build-cmake{,-f,-l} + $SCR/perf-capture.sh baseline .perfeng # full ctest + all bench cases → artifacts/, collated per precision + ``` + `perf-capture.sh` exits non-zero if any precision's baseline isn't fully green → **stop**; + optimization starts only from a clean tree. Build the snapshot table with + `uv run python $SCR/perf-bench.py snapshot artifacts/baseline-bench-

.json --prec

`. + → *Deliverable:* `phase-a-baseline.md` (build config + commit + baseline snapshot table with + a `prec` column) and raw per-precision `artifacts/baseline-{tests,bench}-{d,f,l}.*` — the + exit reference Phase F is judged against. + +- [ ] **[Phase B — pin the correctness net](details/phase-b-correctness-net.md) [HARD GATE].** Inject the smallest fault into + the target (flip an operator, drop a term — agent's call), rebuild, run `build-cmake/tests/checkall`; + `perf-net.py` turns its output into the net (the cases that flip to `-> FAIL`). **No test fails ⇒ the + region is uncovered ⇒ stop and exit** (try a more destructive fault first; then report + the gap — never optimize uncovered code). Revert, confirm green, `git diff` empty. + → *Deliverable:* `phase-b-correctness-net.md` (the fault, the net table + suite + size, + revert confirmed) and `artifacts/fault.diff`; on gate failure this is the `blocked` + coverage-gap report instead. + +- [ ] **[Phase C — pin the performance metric](details/phase-c-performance-metric.md) [HARD GATE].** Measure the target's + benchmark case(s), inject a slowdown (e.g. wrap the body in an N-times loop but beware of compiler optimizations, results still + correct), re-measure. Cases whose `median_ns` rise clearly are your metric. **No + benchmark moves ⇒ the region is uncovered/tooling can't measure ⇒ stop and exit**. Revert, `git diff` empty. + → *Deliverable:* `phase-c-performance-metric.md` (target baseline snapshot, the metric + case(s) with before/after, revert confirmed) and `artifacts/slowdown.diff`; on gate + failure this is the `blocked` no-metric report instead. + +- [ ] **[Phase D — rounding-error analysis](details/phase-d-error-analysis.md).** *Before any speed + edit*, scrutinize the target's rounding error under the **standard model** + (`fl(a∘b)=(a∘b)(1+δ), |δ|≤u`): derive the forward-error bound with any parameter + (`N`, `M`, or other) dependence explicit, and reconcile it with the Phase-A baseline errors (and a + [`scripts/perf-trend.py`](scripts/perf-trend.py) error-vs-`N` fit if growth is in question). + Conclude with a verdict that **sets the inner loop's objective** — `clean` (error at the best + achievable floor → optimize for speed without regressing the derived bound) or `improve-first` + (an *avoidable* `N`/`M` dependence, cancellation, or unstable recurrence → **root out the error + first, then tune speed**; often both improve). Seed the risk table with the size/precision + dependence the analysis exposes. → *Deliverable:* `error-analysis.html` (self-contained, + MathJax — the full derivation) **and** `phase-d-error-analysis.md` (verdict + objective). + +- [ ] **[Phase E — inner loop](details/phase-e-inner-loop.md).** Optimize *toward the Phase-D + objective* against the *narrow* B-net + C-metric only. After every change: rebuild **all three + precision trees** and run the B-net (`checkall`, add `checkall_threads` if OpenMP touched) in + each — it must stay green in float, double, *and* long double; re-measure the C-metric in each — + `median_ns` drops and never rises beyond noise. The net is narrow, so 3× is cheap and catches a + float-only break early. Fast feedback, but not authoritative — Phase F is. Keep the **risk + table** current: for each change, note the accuracy/correctness side effects this narrow net + can't see (size-dependent, precision-specific, input-range — see + [risk-assessment](details/risk-assessment.md)); when a risk is material and cheaply testable, + extend the net ([extending-tests](details/extending-tests.md)) to prove or disprove it. + **Commit each kept change** (and any permanent test addition) as its own self-contained commit + as you go — Phase G squashes them later. + → *Deliverable:* `phase-e-inner-loop.md` (iteration journal: per-change net result + metric + before→after, per precision) and the current `artifacts/change.diff`. + +- [ ] **[Phase F — exit gate](details/phase-f-exit-gate.md).** Re-run the ENTIRE Phase-A baseline + (full `ctest` + ALL benchmark cases) **in all three precisions** — `perf-capture.sh final …` + then `uv run python …/perf-bench.py compare --taskdir …` (deterministic noise-rule verdict). + Exit only when, *in float, double, and long double*: the full suite passes as in Phase A; no + benchmark regresses beyond `max(3·stdev, 2% of median)` *and* the rise survives a re-run; the + target's metric improved or equal; **the Phase-D accuracy objective is honoured** (avoidable + source removed / derived bound not regressed); `git diff` is only the intended change; and the + **risk assessment** is complete and honest — every material risk `proven`, `retired`, + `accepted`, or `residual` ([risk-assessment](details/risk-assessment.md)). A regression in + *any* precision fails the gate. A `proven` material accuracy drop is a **hard no** — fix or + revert, never ship; an unsettled `residual` material drop may land but only as a `partial` run + that surfaces it (with the settling test proposed) for the human reviewer. Any failure ⇒ loop + back to E, or revert and report why. → *Deliverable:* `phase-f-exit-gate.md` (per-precision + comparison table over ALL cases + the six-point checklist + verdict) with raw + `artifacts/final-{tests,bench}-{d,f,l}.*`; then close out the tracker (**Status** = complete | + reverted, **Outcome** one-liner) and write the human report + `summary.html` (the reviewer-facing walkthrough — produced on *every* exit, including a Phase + B/C hard-gate block; it presents every phase's result + numbers, **embeds the required charts** + via `perf-summary.py charts`, and links every deliverable + artifact — `perf-summary.py check` + verifies nothing is orphaned). A faster target bought with a regression or a broken test + elsewhere is **not** a success. + +- [ ] **[Phase G — conclude](details/phase-g-conclude.md).** Hand the work off, in two scripted + steps around the interactive publish. **Squash** the run's intermediate commits into one — + `perf-conclude.sh squash -m ""` (deterministic: preflights a clean tree on a feature + branch, `git reset --soft ` + one commit). Then **offer to push** the branch and + **optionally open a PR** to `develop`, labelled **`perf-eng`** (create the label once if absent) + — both opt-in, confirm with the user. Once the PR number `N` is known, **package**: + `perf-conclude.sh package ` renames `.perfeng/` → `.perfeng-pr-/` (no leading zeros — so a + downloaded archive unpacks into a PR-unique dir) and zips it; **attach that zip** to the PR + (drag-drop in the web UI — `gh` has no binary PR-attach; deliverables are gitignored, so they + never ride in the commit). **If follow-up work lands after conclude** (a reviewer question → + more code + deliverable changes), commit the code, refresh `.perfeng-pr-/`, re-run + `perf-conclude.sh package `, and **overwrite** the zip in its existing PR comment so the + attached archive stays current. → *Deliverable:* one squashed commit; (if published) a + `perf-eng`-labelled PR with the `perfeng-pr-.zip` deliverables archive attached (kept current). + +## Key rules + +- **Measure all three precisions** (float · double · long double) at the baseline, the inner + loop, and the exit gate — one CMake tree per precision (`build-cmake{,-f,-l}`). A + precision-agnostic target can regress in just one precision; double-only misses it. See + [precision-matrix](details/precision-matrix.md). +- **A green suite is narrow, not a guarantee.** Passing tests prove correctness only for + the sizes, grids, and distributions the cases cover. Every run produces a **risk + assessment** of the accuracy side effects that fall outside that coverage — proven or + unproven, all surfaced in `summary.html` — and extends the net to settle the material + ones ([risk-assessment](details/risk-assessment.md), [extending-tests](details/extending-tests.md)). +- **Accuracy is analyzed before speed, not just tested.** Phase D derives the target's + rounding-error bound under the standard model *before any speed edit* and sets the objective: + an `improve-first` verdict makes *removing avoidable error* the inner loop's first job; a + `clean` verdict makes the derived bound — not merely test pass/fail — the accuracy yardstick the + exit gate enforces ([phase-d-error-analysis](details/phase-d-error-analysis.md)). +- **Deterministic steps are scripted, not hand-run.** [`scripts/`](scripts/) drives the build, + capture, JSON collation, test-net parsing, noise-rule comparison, trend fit, chart generation, + and the write-up completeness check so they are reproducible and token-light; the judgement + steps (the fault, the slowdown, the analysis, the optimization) stay with you. See + [scripts/README.md](scripts/README.md). +- **The write-up shows the whole run.** `summary.html` is the reviewer's artifact: every phase's + result with numbers, the **required charts** (performance speedup; error-vs-`N` trend), and + links to *every* deliverable and raw artifact — generated and verified by `perf-summary.py + charts`/`check`. See [deliverables.md](details/deliverables.md#required-visualizations). +- **Walltime is the only metric — local and CI both use it** (CodSpeed Macro Runners run + walltime), so a local result previews the CI gate. No simulation/instruction-count gate, + no CodSpeed account, token, MCP, or cloud upload in the loop. +- **Compare medians, never single runs** (walltime is noisy). Untouched, byte-identical + cases swing several percent run to run — that's noise, not a regression; re-run to confirm. +- **One CMake tree per precision drives that precision** (tests + benchmark), built with + `-DNFFT_BENCHMARK_MODE=walltime`. The `simulation` build is used *only* for the optional + callgrind tie-breaker on an ambiguous untouched control case ([measurement-modes](details/measurement-modes.md)). + Don't use the legacy Autotools `--with-codspeed` path. +- **OpenMP-only changes show only in `checkall_threads`** (it links the `_omp` library). +- Respect the precision-agnostic conventions (`Y()`/`X()`/`FFTW()` mangling, `R`/`C`/`E` + types) — see `CONVENTIONS.md`. Keep the float/double/long-double matrix building. diff --git a/.claude/skills/nfft-perf-eng/details/caveats.md b/.claude/skills/nfft-perf-eng/details/caveats.md new file mode 100644 index 00000000..8245e7e5 --- /dev/null +++ b/.claude/skills/nfft-perf-eng/details/caveats.md @@ -0,0 +1,55 @@ +# Caveats + +*[← Overview & map](../REFERENCE.md) — cross-cutting reference, consult from any phase.* + +- **Crash-class faults lose granularity.** A Phase-B fault that corrupts buffer sizes + or indices (rather than just values) can make a case under-allocate and `SIGABRT` + (exit 134) partway through the suite. The `-> FAIL` lines printed *before* the abort + are still a reliable net, but the XML may be truncated. Prefer a fault that yields + wrong-but-finite results — e.g. the imaginary-sign flip in `trafo_direct`'s kernel + (Phase B) — over one that aborts when you want a complete failure list. +- **OpenMP-only changes show only in `checkall_threads`.** That binary links the + `_omp` library; `checkall` does not. A change guarded by `#ifdef _OPENMP` (or + affecting only parallel scheduling) will leave `checkall` green and must be judged + by `checkall_threads`. +- **A wrong target is possible — but Phase C catches it, so you needn't pre-judge.** + Some code is simply not on any timed path: `init`-only helpers like `intprod()` run + only in `*_init*` (plan setup), and the benchmarks call `init_*` **outside** the + timed `for (auto _ : state)` loop (`benchmarks/bench_nfft_direct.cpp:65`) — they time + only `trafo_direct` / `adjoint_direct`. You don't have to know that up front: the + Phase C hard gate makes it self-evident — the deliberate slowdown moves *no* + benchmark, so there is no metric, so the loop **stops**. The remedy is to target code + inside the timed region (`trafo*` / `adjoint*`) or add a benchmark that covers it. +- **Benchmark coverage is narrow.** Only the **direct** (slow, O(N·M)) transforms + are benchmarked, forward and adjoint, in 1d/2d/3d (`bench_nfft_direct.cpp`). The + **fast** NFFT path (`trafo`, `adjoint`, the `precompute_one_psi` strategies) and + every non-NFFT module have no benchmark yet. Optimizing those requires adding a + benchmark first (Phase C with no coverage is meaningless). +- **`[ERROR] instrument-hooks: failed to write environment.json` is benign.** A + walltime run prints it when the **current working directory** (where the binary + tries to drop `environment.json`) is not writable; the per-case results JSON is + still written to `$CODSPEED_PROFILE_FOLDER/results/` and the exit code is 0. Run + from a writable dir to silence it, or ignore it. +- **Walltime regressions on untouched code are usually noise — but not always random + noise.** A change confined to one function (e.g. `trafo_direct`'s 1d branch) will still + show neighbouring, *byte-identical* cases (`adjoint_direct`, other dims) swinging ±several + percent run to run — worst on the few-iteration 2d/3d cases. This is why Phase F uses a + noise threshold and a re-run, not a bare "nothing got slower" check. +- **A *reproducible* regression on untouched code can still be a code-layout artifact, not a + real slowdown.** Editing the target changes its size, which shifts the alignment of whatever + the linker places after it; an untouched neighbour's hot loop can land on a worse boundary + and read a *stable* +N% that survives re-runs (so the noise-rule re-run does **not** clear + it). The fix is to pin layout: build with `-falign-functions=64 -falign-loops=32` (CI's + `bench-linux.yml` does — the Phase-A flag list now includes them). If a control case still + regresses reproducibly, confirm attribution by re-measuring the **baseline** under the same + aligned flags (`git stash` the change, rebuild, measure): a layout artifact vanishes; a true + coupling persists. Observed concretely: a float `adjoint_direct_1d[32/100]` control read +7.5% + without the flags and −0.9% (flat) with them. +- **A green suite is narrow, not a guarantee.** The Phase-B net is a deliberate subset, + and even the full Phase-A suite only checks specific transforms at specific sizes and + grids. An optimization can be faster, pass every test, and still have dropped accuracy + on a path no case exercises — most often at transform sizes larger than any tested + case, or for node/coefficient distributions no case spans. Passing tests is + authoritative *for what the tests cover*; the gap is what the **risk assessment** + ([risk-assessment.md](risk-assessment.md)) surveys and [extending-tests.md](extending-tests.md) + can close. diff --git a/.claude/skills/nfft-perf-eng/details/deliverables.md b/.claude/skills/nfft-perf-eng/details/deliverables.md new file mode 100644 index 00000000..275742fd --- /dev/null +++ b/.claude/skills/nfft-perf-eng/details/deliverables.md @@ -0,0 +1,246 @@ +# Deliverables & the `.perfeng/` directory + +*[← Overview & map](../REFERENCE.md) — cross-cutting reference; read before Phase A, consult from every phase.* + +This loop is **front-to-back auditable**: every phase produces a concrete +deliverable, and a phase is not exitable until its deliverable exists in the +directory and the tracker row is updated. The tracker is the live front page; the phase +docs are the durable record of what was measured and decided. + +**Deliverables are never committed.** They live in a fixed, **gitignored** directory, +`.perfeng/`. Git only ever sees your *source* changes — so commit each self-contained unit of +work as you go (a kept optimization step, a permanent test addition), and let +[Phase G](phase-g-conclude.md) squash them into one commit and carry the deliverables to the +reviewer as a **zip attached to the PR**. The run's permanent record is that commit + PR + +archive, not files in the tree. + +**The rule:** *deliverable = exit gate.* You may not mark a phase `done` (or move to +the next phase) until its deliverable file is written and the tracker reflects it. +For the hard-gate phases (B, C) the deliverable also records the gate verdict — and +if the gate fails, the deliverable is the **stop report** (status `blocked`), not a +green light. + +**Don't hand-write the deliverables — fill the templates.** Every deliverable below +has a fill-in skeleton under [`../templates/`](../templates/) (one per phase, plus the +tracker). At each phase, copy the matching template into `.perfeng/` and fill its +`<…>` placeholders; the canonical tables below are pre-embedded in them. This is what +keeps runs comparable. The phase docs name the exact template to copy. + +## Where it lives + +One fixed, **gitignored** directory — `.perfeng/` — holds the **current** run. It is a +*static* path (no `NNNN` sequence, no shared index): each run overwrites it, and Step 0's +[`scripts/perf-init.sh`](../scripts/perf-init.sh) moves any prior `.perfeng/` aside to +`.perfeng.bak/` (a one-level undo) before recreating it. Nothing here is committed. At +[conclude](phase-g-conclude.md), once the PR number `N` is known, `.perfeng/` is renamed to +`.perfeng-pr-/` (no leading zeros) and zipped, so a downloaded results archive unpacks into a +PR-unique directory. (`.perfeng-pr-*` is gitignored too.) + +``` +.perfeng/ # gitignored; the current optimization, front-to-back + BASE # squash base — the commit + branch HEAD was at, at Step 0 (for Phase G) + README.md # THE TRACKER — status per phase, links, outcome (front page) + phase-a-baseline.md # Phase A deliverable — baseline snapshot (the exit reference) + phase-b-correctness-net.md # Phase B deliverable — the pinned net (or blocked report) + phase-c-performance-metric.md # Phase C deliverable — the pinned metric (or blocked report) + phase-d-error-analysis.md # Phase D deliverable — rounding-error verdict & objective (companion to the HTML) + error-analysis.html # Phase D deliverable — full standard-model analysis (MathJax, self-contained) + phase-e-inner-loop.md # Phase E deliverable — iteration journal + phase-f-exit-gate.md # Phase F deliverable — final verdict + summary.html # Close-out — human-facing report (any outcome) + artifacts/ # raw captured data, kept verbatim for exact diffing + baseline-tests-{d,f,l}.log # Phase A (tee'd ctest output, one per precision) + baseline-bench-{d,f,l}.json # Phase A (per-case stats, collated by perf-capture.sh) + fault.diff # Phase B (the injected correctness fault) + fault.log # Phase B (checkall stdout under the fault — parsed by perf-net.py) + net-cases.txt # Phase B (failing-case names = the net; revert + inner-loop checks) + slowdown.diff # Phase C (the injected slowdown) + trend-{baseline,optimized}-{d,f,l}.dat # Phase D (differential trend study: N→error pairs, per precision, if run) + change.diff # Phase E/F (the actual optimization) + final-tests-{d,f,l}.log # Phase F + final-bench-{d,f,l}.json # Phase F (collated like baseline-bench-*.json) + chart-speedup-{d,f,l}.svg # Phase F (perf-summary.py — embedded tabbed in summary.html) + chart-trend-{d,f,l}.svg # Phase F (perf-summary.py — when a trend study ran, per precision) +``` + +Captures are **per precision** — `-d` (double), `-f` (float), `-l` (long double) — because +Phases A, D, E run the whole float·double·long-double matrix; see +[precision-matrix.md](precision-matrix.md). (`fault.diff`/`slowdown.diff`/`change.diff` are +source edits, precision-independent, so they are single files.) + +**Test additions.** When the [risk assessment](risk-assessment.md) leads you to extend +the net ([extending-tests.md](extending-tests.md)), a *permanent* case lands as ordinary +source: harness/registration edits under `tests/` and any committed `tests/data/*.txt` + +generated headers — part of the change's `git diff` (commit it), not under `.perfeng/artifacts/`. +A *temporary* probe is removed before close-out; only its finding survives, in the risk table. + +The benchmark binary writes one file per process to +`$CODSPEED_PROFILE_FOLDER/results/.json` (a transient scratch dir). +[`scripts/perf-capture.sh`](../scripts/perf-capture.sh) collates these into the **kept** +deliverable — one flat-array JSON per precision under `artifacts/` +(`baseline-bench-{d,f,l}.json` in Phase A, `final-bench-{d,f,l}.json` in Phase F). +[`scripts/perf-bench.py compare`](../scripts/perf-bench.py) diffs the two, per precision, +applying the noise rule. + +- The phase docs mirror the skill's own `details/` phase names 1:1, so the mapping + from instruction to deliverable is obvious. +- Keep `artifacts/` raw and verbatim (logs, JSON, diffs) — the narrative docs embed + *summaries* (tables) of these; the raw files are what Phase F diffs against. +- Redirect commands that the phase docs write to `/tmp/...` into `.perfeng/artifacts/` + instead, so the capture is durable. Example: `... | tee + .perfeng/artifacts/baseline-tests.log`. + +## Step 0 — open the `.perfeng/` directory (gated) + +Bootstrap the directory first. This is deterministic — run +[`scripts/perf-init.sh `](../scripts/perf-init.sh): it creates the fixed +`.perfeng/` (with `artifacts/`), copies **every** template in (tracker→`README.md`, the +phase docs, and `error-analysis.html`), stamps the baseline commit, and records the **squash +base** — the commit + branch HEAD points at right now — into `.perfeng/BASE` so +[Phase G](phase-g-conclude.md) can collapse the run's commits deterministically. If `.perfeng/` +already exists it is moved aside to `.perfeng.bak/` and recreated fresh. What remains manual: + +1. **Set the tracker's Target line** — open `.perfeng/README.md` (copied from + [`../templates/tracker.md`](../templates/tracker.md)) and set the **Target** line; leave + **Outcome** as `—` and **Status** as `in-progress`. + +*Deliverable = exit gate:* Step 0 is not done until `.perfeng/README.md` exists with all six +phase rows present and `.perfeng/BASE` records the squash base. Every later "flip the tracker +row" instruction assumes this infrastructure already exists. + +## The tracker (`README.md`) + +The front page of the directory. Created in Step 0 by copying +[`../templates/tracker.md`](../templates/tracker.md), updated at every phase boundary: +flip the row's status, fill its **Deliverable** link and **Exit signal**, and append a +**Log** line. The template ships with **Phase A** `🔄` (everything below `⬜`, `Outcome` +undecided) so a fresh copy already reflects Step 0's end state — Phase A is the first phase. + +Status legend: `⬜` todo · `🔄` in-progress · `✅` done · `⛔` blocked (hard gate +failed — see that phase's deliverable for why). **Lifecycle:** a row is `⬜` until +entered, `🔄` on entry, `✅` on exit-gate pass, `⛔` on a hard-gate stop or abandon. +**Header mapping:** a `⛔` on any row ⇒ header **Status** = `reverted`. When the loop +ends, set the header **Status** to `complete` (Phase F passed) or `reverted` (gave up +/ hard gate) and fill **Outcome**. **Reopen:** if Phase F bounces the work back +([phase-f-exit-gate.md](phase-f-exit-gate.md)), reset the Phase E row from `✅` to `🔄` and add a +fresh journal row — `✅` means "passed *its own* gate", not "final". + +## The human summary (`summary.html`) + +The markdown deliverables are the audit trail for the next agent; `summary.html` is the +report a **human reviewer** reads to approve the work. Produce it from +[`../templates/summary.html`](../templates/summary.html) at **close-out — on every +exit**: a completed run (Phase F passed), a reverted run (gave up at Phase F), or a +hard-gate block (Phase B or C). It walks the reviewer concisely from the idea (what to +change and why), through the process, to the results — and on a failure documents the +known and hypothetical causes. It also banks **patterns used/discovered** and +**carry-forward notes** so knowledge accumulates across optimization targets. Set its +`` to `ok` / `partial` / `fail` to match the outcome. In [Phase +G](phase-g-conclude.md) `summary.html` (with the rest of `.perfeng/`) is zipped and attached +to the PR — it is the reviewer's entry point into the evidence bundle. + +The run's *other* human-facing HTML is the Phase-D **`error-analysis.html`** — the +self-contained, MathJax-rendered rounding-error analysis (the math the reviewer must follow +to trust the accuracy verdict). It is produced mid-run (before any speed edit), not at +close-out, and `summary.html` links to it rather than repeating the derivation. + +`summary.html` must present the **entire** run to the reviewer — every phase's result, with +the numbers, **and links to every deliverable and raw artifact** (the per-phase docs, +`error-analysis.html`, and each file under `artifacts/`: tests logs, benchmark JSON, the +diffs, the net, the charts). Nothing in the directory may be orphaned; +[`scripts/perf-summary.py check`](../scripts/perf-summary.py) enforces this at close-out +(exit 1 on any unlinked file or missing required chart). + +### Required visualizations (the lower bar) + +`summary.html` must *show*, not only tabulate, the results a reviewer judges the run on. +Generate these deterministically with +[`scripts/perf-summary.py charts`](../scripts/perf-summary.py) — self-contained inline SVG, +no JS or network — and embed them (``). This is the +**minimum**; richer/extra charts are welcome. + +1. **Performance — % faster per case, one chart per precision** + (`artifacts/chart-speedup-{d,f,l}.svg`, shown **tabbed** in `summary.html` so the full case + list stays readable). **Required on every completed run.** Base→final per case; green + faster · red regressed · grey within-noise, with the noise threshold drawn. Tab across + float·double·long double to see whether the win holds everywhere and whether any control + case regressed. (A correctness-only precision has no chart — drop its tab.) +2. **Accuracy trend — error vs `N`, log-log, baseline vs optimized, one chart per precision** + (`artifacts/chart-trend-{d,f,l}.svg`, **tabbed** in `summary.html`). When a differential + trend study is run at all (any `improve-first` verdict, or a size-dependent accuracy risk), + **run it in all three precisions** — float and long double can behave differently (float + cancels/overflows sooner), so a double-only trend would hide a precision-specific accuracy + change. Each chart plots both fitted exponents `p` (expect ≈0.5, √N) side by side — the + clearest single picture of an order-of-growth change. +3. **Accuracy margin — measured error vs the documented bound**, per case × precision (table + always; a bar chart optional). Show the headroom to `C·ε`; if the change touched accuracy, + show error before→after. + +On a **Phase B/C hard-gate block** there is no final benchmark/trend data, so the charts are +`N/A` — `summary.html` then carries the §"Why it stalled" analysis instead, and +`perf-summary.py check` does not demand the charts. + +## Canonical formats + +Use these exact shapes so snapshots are comparable across phases (Phase F diffs the +Phase-A and final benchmark snapshots; the B net is re-checked in E and F). + +**Benchmark snapshot** — Phase A, C, E, F. Embed a table in the narrative doc; keep +the raw per-case JSON in `artifacts/`. The `prec` column (`d`/`f`/`l`) keeps all three +precisions in one comparable table ([precision-matrix](precision-matrix.md)). + +```markdown +| prec | case | median_ns | stdev_ns | rounds | +|------|----------------------------|-----------|----------|--------| +| d | nfft_forward_direct_1d/… | 123456 | 789 | 50 | +| f | nfft_forward_direct_1d/… | 98765 | 654 | 50 | +| l | nfft_forward_direct_1d/… | 210987 | 912 | 50 | +``` + +**Correctness net** — Phase B. The set of cases that flip to `-> FAIL` under the +injected fault, plus the suite to run in the inner loop and the net size (the latter +two as header bullets above the table, as in the template): + +```markdown +- **Suite to run in inner loop:** `nfft` +- **Net size:** 149 cases + +| suite | case | error | bound | +|-------|----------------------------------------|---------|----------| +| nfft | nfft_1d_50_50.txt … trafo_direct | 5.7e+14 | 1.07e-14 | +``` + +**Comparison table** — Phase F. Baseline vs final per case **per precision**, with the +noise rule applied and a per-case verdict. A regression in *any* precision fails the gate. + +```markdown +| prec | case | base median_ns | final median_ns | Δ% | threshold | verdict | +|------|--------------------------|-----------------|-----------------|------|-----------|---------| +| d | nfft_forward_direct_1d/… | 123456 | 95012 | −23% | 2%/3σ | ✅ faster | +| f | nfft_forward_direct_1d/… | 98765 | 80120 | −19% | 2%/3σ | ✅ faster | +| l | nfft_forward_direct_1d/… | 210987 | 165430 | −22% | 2%/3σ | ✅ faster | +``` + +**Iteration journal** — Phase E. One row per change attempt, appended as you go; +`net` is `green` or the B-net case that flipped to `-> FAIL`; the median columns are +the Phase-C metric case(s), before→after (medians, never single runs). + +```markdown +| iter | change | net | metric median (ns) before→after | +|------|--------------------------------|----------------------------------|----------------------------------| +| 1 | hoist `K[j]` out of inner loop | green | 123456 → 118900 | +| 2 | unroll ×4 | FAIL: nfft_1d_50_50 trafo_direct | — (reverted) | +``` + +**Risk table** — Phase B (seed), D (analysis seeds), E (per-iteration rows), F (consolidated). The negative +side effects the narrow net can't see, each with its category, state, and disposition. +States: `proven` (fix/revert) · `retired` (a check that would expose it now passes) · +`accepted` (within the documented error bound / non-accuracy only — never a material drop) · +`residual` (suspected material drop, no constructible check; surfaced in `summary.html`). +See [risk-assessment.md](risk-assessment.md). + +| risk | category | state | evidence / disposition | +|---------------------------------------------|-------------------|----------|-------------------------------------------------| +| accuracy may drop for N≫ tested 1d sizes | size-dependent | retired | added online check N=4096, green d/f/l (temp probe) | +| reassociated sum loses bits in long double | accuracy-for-speed| residual | plausible; no case isolates it — see summary | +| shared `K[j]` helper reused by `nnfft` | out-of-scope coupling | accepted | `nnfft` suite green in Phase F; within bound | diff --git a/.claude/skills/nfft-perf-eng/details/extending-tests.md b/.claude/skills/nfft-perf-eng/details/extending-tests.md new file mode 100644 index 00000000..eed7b3b7 --- /dev/null +++ b/.claude/skills/nfft-perf-eng/details/extending-tests.md @@ -0,0 +1,146 @@ +# Extending the net (proving or disproving a risk) + +*[← Overview & map](../REFERENCE.md) — cross-cutting reference; the response to a risk +([risk-assessment.md](risk-assessment.md)) the existing net can't settle.* + +When the [risk assessment](risk-assessment.md) surfaces a **material** risk the current +net can't see — most often *"accuracy may have dropped at sizes larger than any tested +case"* — the right move is often to **add a test that would expose it**, run it, and +convert the risk from `residual` to `retired` (or to `proven`, which means fix/revert). +This is a judgement call, not a mandatory gate: extend the net when the risk is plausible +*and* matters *and* a probing check is cheap relative to the doubt it removes. + +This doc is the *how*. The test harness and the two test classes it builds on are +documented in [`test-methodology.md`](../../../docs/agents/test-methodology.md) — read it +first; this doc only adds the perf-loop-specific framing. + +## The two levers (glossary class names) + +`CONTEXT.md` fixes the two CUnit test classes as the **file-based check** and the +**online check** (see also `test-methodology.md`). They give two ways to extend coverage, +with very different costs: + +1. **Online check (fast-vs-direct) — cheap, no data file.** `setup_online` / + `setup_adjoint_online` generate random input at a chosen size, build the reference by + running the C **direct** transform, and compare the **fast** transform against it. To + add one you register another case in the test harness (a new size in the online list + for the suite) — **no `tests/data/*.txt` file is needed**. This is the primary lever + for settling a **fast-path** optimization's risk at a larger `N`/`M`: the trusted + direct transform *is* the oracle — **but only when that direct transform is itself + sound** (its own error already at the floor, per the [Phase-D analysis](phase-d-error-analysis.md)). + If you are optimizing the direct transform, it cannot also be its own oracle — use lever 2. +2. **File-based check (reference-data generator) — adds a high-precision oracle.** A + `tests/data/*.txt` reference is produced offline by the **reference-data generator** + (`tests/refgen/`, "refgen"), and both the direct and fast transforms are checked + against it. Use this when you need an *independent* oracle — chiefly when the **direct** + transform itself is the target (an online check would compare the direct transform + against itself and prove nothing), or for a permanent regression pin at a specific + grid. Generate with: + + ```bash + uv run --with mpmath==1.3.0 python -m tests.refgen.generate --module --precision 64 + ``` + + (NFFT/NFCT/NFST only; see `test-methodology.md` for grids, options, and the + committed-artifact workflow. Per [[prefer-uv-not-pip]], always `uv run` — never pip.) + +## Which lever — by what you're optimizing + +The lever differs by optimization type because the *oracle* differs: a fast-path change +can be pinned by the direct transform (an online check), but a direct-path change needs an +external (refgen) oracle. + +| Target of the optimization | Cheapest settling lever | +|----------------------------|-------------------------| +| **fast path** (`trafo`, `adjoint`, `precompute_one_psi`) | online check at a larger size — direct is the oracle, no data file | +| **direct transform** (the worked example) | file-based check at a new grid (refgen) — an online check would self-compare | +| **OpenMP-only** change | the same case in `checkall_threads` (links the `_omp` lib); a serial-only case won't see it | +| **non-NFFT module / no benchmark** | the module's own suite; coverage and even the benchmark may be missing first (see [caveats.md](caveats.md)) | + +## Temporary vs permanent + +Decide this explicitly and record it in the risk table: + +- **Temporary probe** — added solely to settle one risk for this run; **prefer an online + check** (no data files to manage). Run it on the *baseline* (unoptimized) code to + confirm it passes, then on the change; the delta is the evidence. **Remove it before + close-out** so Phase F's "`git diff` is only the intended optimization" exit condition + holds — and record in the risk table that it was a temporary probe and what it showed. +- **Permanent check** — kept in the suite and CI; this is usually where a *file-based* + (refgen) addition lands, since generated data is awkward to add and remove for one run. + Justified when the gap it closes is real and its run cost is acceptable (online checks at + moderate sizes are cheap; large refgen grids add to every CI run — weigh it). A permanent + addition is a *legitimate part of the change's `git diff`*: harness/registration edits, + and any committed `tests/data/*.txt` + generated headers. **Confirm it green on the + baseline commit** (it pins existing behaviour at a new condition, so it must pass + unoptimized) — that makes it a known-green case Phase F compares against alongside the + Phase-A baseline. Note it in the Phase-F deliverable so the reviewer expects test files. + +## The settle loop (proving or disproving) + +1. Build the settling check (online or refgen) at the size/distribution the risk names. +2. Confirm it **passes on the baseline** code — `git stash` the change (or check out the + baseline commit), build, run the new case. If it fails on baseline, the check is + mis-specified (or the bound is wrong), not evidence about your change — fix the check. +3. Restore the change, run the case again **in all three precisions**: + - passes → risk **retired** (record it; remove if temporary). + - fails → risk **proven** → for a material accuracy drop this is a **hard no**: fix the + optimization or revert (Phase E / F). Never ship it. +4. Record the outcome in the risk table ([deliverables.md](deliverables.md#canonical-formats)). + If no settling check is constructible after honest attempts, the risk stays `residual` + → the run lands as `partial` with the gap surfaced (see [risk-assessment.md](risk-assessment.md)). + +## Differential trend analysis (the strongest settle for accuracy-for-speed) + +A single pass/fail at one size only shows the change is *within the bound there* — it cannot +tell you whether the optimization changed the error's **order of growth**. An optimization can +sit under the bound at every committed size yet still be heading for trouble at larger N (e.g. +a reassociation or recurrence that turns O(√N) round-off into O(N)). The check the existing +tests can't give: **the trend with respect to a parameter, compared against the unmodified code.** + +So for a size-dependent / accuracy-for-speed risk, prefer to *also* run a **differential trend +study**: + +1. Pick a geometric sweep of the parameter (e.g. `N = 32, 64, …, 4096`), with **one + higher-precision oracle per point** (the [file-based check](#the-two-levers-glossary-class-names) + — a direct-transform target needs an *external* oracle, since it cannot validate itself). + > **The oracle must be strictly more precise than the build under test — use refgen + > (mpmath, arbitrary `dps`), never a same-precision C reference.** Computing the + > reference in C `long double` works for the *double* and *float* builds (quad > double > + > float), but it is **not an oracle for the long-double build itself** — a quad reference + > cannot out-resolve a quad transform, so the long-double error it reports is meaningless + > (it measures the recurrence's *added* error vs an exact-trig long-double sum, not the + > true error). Drive `tests.refgen.transforms` directly at the swept `N`/`M` with + > `mpmath.mp.dps` set comfortably above the build's decimal digits (e.g. `dps ≥ 50` covers + > double; `dps ≥ 40` covers quad) — the *same* arbitrary-precision oracle serves all three + > precisions uniformly. mpmath is slow, so shrink `M` or thin the sweep if a point is + > expensive; that is cheaper than maintaining a parallel C oracle and is correct in every + > precision. +2. Measure the suite's error metric at each point for **both** the optimized code and the + **baseline** (`git stash` the change between runs), **in each precision** (float and long + double can diverge — float cancels/overflows sooner, so a double-only sweep hides a + precision-specific change). Keep the other parameters fixed (e.g. `M`) so only the swept one + varies. (The harness `MAX_SECONDS` cost guard skips large direct cases — raise it temporarily + for the sweep, or shrink `M`; revert after.) +3. Fit `error ~ C·Nᵖ` (log-log least-squares) for each precision with + [`scripts/perf-trend.py`](../scripts/perf-trend.py) and **compare the exponents** per + precision. **The verdict is Δp (optimized − baseline), not the absolute `p`.** The absolute + exponent often sits *above* the textbook `√N` (`p≈0.5`): a sub-dominant working-precision + term — e.g. forming the phase `2π(k−N/2)x_j` carries absolute error `~N·u` — lifts the + measured `p` to ≈0.6–0.7 even for the *baseline* code, which is fine and expected. What + matters is that the optimization does not *steepen* it: a Δslope near zero (both arms on the + same curve) ⇒ order of growth preserved ⇒ risk **retired**, with far more confidence than any + single-size pass. A clearly steeper optimized slope (in *any* precision) ⇒ a **proven** + order-of-growth regression (hard no), even if every individual point is still under the flat + bound. +4. The sweep grids/oracles are a *temporary probe* (revert the grid/cost edits at close-out); + keep the per-precision `(N, error)` data as `artifacts/trend-{baseline,optimized}-{d,f,l}.dat` + ([`perf-summary.py charts`](../scripts/perf-summary.py) charts them tabbed in `summary.html`) + and cite the exponents per precision in the risk table. (Worked example: the `trafo-direct` + reference run — blocked + recurrence vs baseline: Δp = +0.034 (double), +0.0012 (float), +0.005 (long double), all + ≪ the ±0.15 tolerance ⇒ retired. Absolute `p` ≈ 0.62–0.68, above √N, in *both* arms — the + working-precision phase term, not a recurrence artifact; the small Δp is what settles it.) + +See [`test-methodology.md`](../../../docs/agents/test-methodology.md) for the harness +mechanics (case registration, `NN` formula, bound tables, adding a whole new transform). diff --git a/.claude/skills/nfft-perf-eng/details/measurement-modes.md b/.claude/skills/nfft-perf-eng/details/measurement-modes.md new file mode 100644 index 00000000..75a03c9b --- /dev/null +++ b/.claude/skills/nfft-perf-eng/details/measurement-modes.md @@ -0,0 +1,72 @@ +# Measurement: walltime everywhere (local and CI) + +*[← Overview & map](../REFERENCE.md) — cross-cutting reference, consult from any phase.* + +**There is one performance metric: wall-clock time.** CI measures it with **CodSpeed +Macro Runners** in walltime mode; you measure it locally with the same `walltime` build, +run directly and offline. Local and CI agree on *what* is measured, so a local result is +a faithful (if noisier) preview of the CI gate. There is no separate simulation / +instruction-count gate to satisfy, and no CodSpeed account, token, MCP, or cloud upload +in the loop. + +Mode is **baked in at build time** via `-DNFFT_BENCHMARK_MODE=` (`off`, the default, +builds no benchmarks): + +| | **walltime** — the metric, local and CI | +|---|---| +| Build | `-DNFFT_BENCHMARK_MODE=walltime` | +| Run | the binary directly, offline (CI runs it on a CodSpeed Macro Runner) | +| Metric | wall-clock ns: min/mean/median/stdev/IQR | +| Local result | JSON at `$CODSPEED_PROFILE_FOLDER/results/.json` | +| Noise | real timing noise — compare **medians**, re-run to confirm | + +The benchmark build still uses codspeed-cpp (FetchContent fetches it; submodules +auto-recurse) — that is just the harness library and needs **no account**. The whole +loop, Phases A/C/E/F, reads the local walltime JSON. + +## The noise rule (this is the metric) + +Because walltime is noisy in a container, **manage the noise explicitly** — this is the +metric, so the discipline below *is* the gate, not a workaround: + +- Compare **medians**, never single runs. +- Treat a case as regressed only when its median rises **beyond noise** — past, say, + `max(3·stdev, 2% of the median)`. Identical, untouched code routinely swings several + percent here (worst on the few-iteration 2d/3d cases), so a hard "no case may be + slower at all" rule produces **false regressions**. +- **Re-run** any case that trips the threshold before believing it; noise rarely + survives a second run, a real regression does. This is scripted: + [`perf-confirm.sh `](../scripts/perf-confirm.sh) re-measures only the flagged + precision(s) and reports which regressions survive. Re-measure on a **quiet** machine — + `perf-capture.sh` logs the 1-min loadavg and WARNs when it is high, because concurrent + builds/sweeps inflate medians (a whole capture can read +10% under load yet sit within + noise when re-run idle). Use `perf-capture.sh … --bench-only` to re-measure without the + slow deterministic `ctest`. + +## Deterministic tie-breaker (optional) — is a stuck control case real work or layout? + +When an **untouched** control case stays flagged even after a quiet re-run, the question +is *why*: real coupling to your change (extra work executed) or a **code-layout/cache +artifact** (your edit shifted the alignment of an untouched neighbour, so the same +instructions now run on a worse cache/branch boundary). Wall-clock alone cannot tell +these apart; the deterministic instruction count can: + +- Build a one-off `simulation` tree (`perf-build.sh simulation`, or just one precision) + and run the affected case under `valgrind --tool=callgrind … + --benchmark_filter=''`. Read `I refs` from the callgrind output. (This is a + process **total** that includes fixed startup ~3M Ir, so use it as a *before/after + delta* on a single case, not an absolute.) +- **Same `I refs` before and after, but a higher wall-clock ⇒ layout/cache, not real + work.** The fix is to **pin layout** with the CI-matching alignment flags + (`-falign-functions=64 -falign-loops=32`, already in the Phase-A flag list and in + `bench-linux.yml`) so the artifact does not appear — *in your run or in CI*. Note: CI + is **also** walltime now, so a layout artifact is **not** invisible to the gate the way + it would have been under an instruction-count gate — it must be mitigated (aligned + flags) or attributed, not waved through. See [caveats.md](caveats.md) (layout regressions). +- **`I refs` rose ⇒ real extra work** — a genuine regression to attribute and fix, + regardless of layout. + +This callgrind cross-check is the **only** use of the `simulation` build mode in the +loop — a diagnostic for an ambiguous control case, not a routine metric and not a gate. +Clean per-case instruction counts (`codspeed run`) need a CodSpeed token and are **not** +part of this loop; the raw single-case `I refs` delta above needs no account. diff --git a/.claude/skills/nfft-perf-eng/details/phase-a-baseline.md b/.claude/skills/nfft-perf-eng/details/phase-a-baseline.md new file mode 100644 index 00000000..d9a369a6 --- /dev/null +++ b/.claude/skills/nfft-perf-eng/details/phase-a-baseline.md @@ -0,0 +1,92 @@ +# Phase A — build tree + full baseline (the exit reference) + +*[← Overview & map](../REFERENCE.md) · Prev: [Step 0 — .perfeng directory](deliverables.md#step-0--open-the-perfeng-directory-gated) · Next: [Phase B — correctness net](phase-b-correctness-net.md)* + +**One self-contained CMake tree per precision drives the loop** — float, double, and long +double (the sources are precision-agnostic, so a change can pass in one precision and break +another; see [precision-matrix](precision-matrix.md)). Each tree is built in **walltime** +mode (`-DNFFT_BENCHMARK_MODE=walltime`) — the one metric, the same mode CI runs on its +CodSpeed Macro Runners ([measurement-modes](measurement-modes.md)). Each tree produces the +library, the CUnit tests, **and** the benchmark binaries. (Don't use the Autotools +`--with-codspeed` path — it is legacy.) This build + the capture below is also what +**confirms the walltime harness is operable** end to end (codspeed-cpp fetched, the binary +writes its per-case JSON) — no separate preflight; if the harness can't run, it fails here. + +The build is **deterministic**, so it is driven by [`scripts/perf-build.sh`](../scripts/perf-build.sh) +— it clears a stale in-source `config.h` (a leftover Autotools double config otherwise +poisons the float/long-double trees, mis-linking `nfft_*` vs `nfftf_*`; CI does the same), +then configures + builds `build-cmake{,-f,-l}` with the **CI-aligned flags** +(`-O3 … -falign-functions=64 -falign-loops=32` — the alignment pair must match +`.github/workflows/bench-linux.yml` so an untouched neighbour can't show a phantom alignment +regression; see [caveats](caveats.md)): + +```bash +SCR=.claude/skills/nfft-perf-eng/scripts +$SCR/perf-build.sh walltime # configures + builds all three precision trees +# configures + builds build-cmake (double), build-cmake-f (float), build-cmake-l (long double) +``` + +(The flags, the `config.h` clearing, and the per-precision `-DNFFT_ENABLE_{FLOAT,LONG_DOUBLE}` +toggles live *inside* the script, the single source of truth — read it for the exact recipe. +Reuse one fetched codspeed-cpp checkout across trees with `PERF_CODSPEED_SRC=`.) + +Both signals come from each tree: + +- **correctness** — `ctest --test-dir `, or run `/tests/checkall` directly for the + granular `-> OK/FAIL` stdout (same CUnit suite as `make check`, at precision-appropriate + tolerances — the very thing that catches a precision-specific break). +- **performance** — `/benchmarks/bench_nfft_direct[_omp]`, run directly (walltime) — the + same metric CI gates on ([Measurement modes](measurement-modes.md)); compare medians under the + noise rule. + +(Autotools `make check` is a valid, CI-canonical correctness path — see +[`test-methodology.md`](../../../../docs/agents/test-methodology.md) — but it does not build the +benchmarks, so the loop stays in the CMake trees.) + +Now record the complete state of **both** signals on each unmodified tree — the contract the +finished work is judged against in [Phase F](phase-f-exit-gate.md). It must be the *full* suite, +not the scoped subset used later, **in every precision**. The capture is deterministic too — +[`scripts/perf-capture.sh`](../scripts/perf-capture.sh) runs the full `ctest` and the full +benchmark per precision straight into the task dir's `artifacts/`, collating the per-process +codspeed JSON into one flat array per precision: + +```bash +$SCR/perf-capture.sh baseline .perfeng +# → artifacts/baseline-tests-{d,f,l}.log and artifacts/baseline-bench-{d,f,l}.json +# exit 0 = every precision fully green and captured; exit 1 = a precision was not (see WARNs) +``` + +If **any** precision's baseline is not fully green (`perf-capture.sh` exits non-zero), **stop** — +optimization starts only from a clean tree. Keep these artifacts for the whole task. + +## Deliverables (exit criteria) + +`phase-a-baseline.md` is **the exit reference** — [Phase F](phase-f-exit-gate.md) is +judged against it, so it must be complete and durable (the `artifacts/` captures, not +vanishing `/tmp`). `perf-capture.sh` already wrote the durable artifacts; build the snapshot +table from them with [`scripts/perf-bench.py snapshot`](../scripts/perf-bench.py): + +```bash +for p in d f l; do uv run python $SCR/perf-bench.py snapshot \ + .perfeng/artifacts/baseline-bench-$p.json --prec $p; done +``` + +Fill [`../templates/phase-a-baseline.md`](../templates/phase-a-baseline.md), which +records: the build config (the `perf-build.sh` walltime mode + the CI-aligned flags), the +baseline commit SHA, the full-suite test result summary **for each precision**, +and the baseline benchmark snapshot as a **Benchmark snapshot** table (canonical format with +the `prec` column — [deliverables.md](deliverables.md#canonical-formats)) covering **all** cases +in **all three precisions** (the `perf-bench.py snapshot` rows drop straight in). + +`artifacts/` (verbatim, per precision): `baseline-tests-{d,f,l}.log` (tee'd `ctest`) and +`baseline-bench-{d,f,l}.json` (per-case stats, collated by the script). The narrative doc embeds +the *summary* table; these raw files are what Phase F diffs against. (If a precision can't +build/run the **benchmark** in your environment, `perf-capture.sh` captures its tests only and +WARNs — record that, never silently drop a precision; see [precision-matrix](precision-matrix.md).) + +**Exit gate** — *deliverable = exit gate*: Phase A is not exitable until the baseline is +**fully green in float, double, and long double**, `phase-a-baseline.md` + the per-precision +raw artifacts exist, and the tracker Phase A row reads `✅` with exit signal `full suite green +(d/f/l); N cases captured`. A non-green baseline in *any* precision ⇒ stop; do not proceed to Phase B. + +*[← Overview & map](../REFERENCE.md) · Prev: [Step 0 — .perfeng directory](deliverables.md#step-0--open-the-perfeng-directory-gated) · Next: [Phase B — correctness net](phase-b-correctness-net.md)* diff --git a/.claude/skills/nfft-perf-eng/details/phase-b-correctness-net.md b/.claude/skills/nfft-perf-eng/details/phase-b-correctness-net.md new file mode 100644 index 00000000..92de7058 --- /dev/null +++ b/.claude/skills/nfft-perf-eng/details/phase-b-correctness-net.md @@ -0,0 +1,80 @@ +# Phase B — pin the correctness net + +*[← Overview & map](../REFERENCE.md) · Prev: [Phase A — baseline](phase-a-baseline.md) · Next: [Phase C — performance metric](phase-c-performance-metric.md)* + +**A1. Identify the target.** A specific function / region, e.g. `X(trafo_direct)()` — +the direct, O(N·M) NDFT in `kernel/nfft/nfft.c:145`. + +**A2. Inject a fault.** Make the *smallest* edit that changes the target's behaviour +— flip an operator, drop a term. The goal is to make dependent tests fail, not to be +realistic. Example: in `trafo_direct`'s 1d branch, flip the sign of the imaginary +kernel — `v += f_hat[k_L] * (COS(omega) - II * SIN(omega));` → +`v += f_hat[k_L] * (COS(omega) + II * SIN(omega));`. + +**A3. Rebuild and see what flips to FAIL.** This set is your correctness net. The fault was +yours to design (it depends on the code); turning the resulting `checkall` stdout *into the +net* is deterministic, so [`scripts/perf-net.py`](../scripts/perf-net.py) does it — the +canonical table and the failing-case name set, no hand-grepping: + +```bash +SCR=.claude/skills/nfft-perf-eng/scripts +D=.perfeng/artifacts +cmake --build build-cmake -j >/dev/null # rebuild the library +build-cmake/tests/checkall > $D/fault.log 2>&1; echo "exit=$?" +uv run python $SCR/perf-net.py table $D/fault.log # the Correctness-net table (suite|case|error|bound) +uv run python $SCR/perf-net.py names $D/fault.log > $D/net-cases.txt # the failing-case set (revert + inner-loop checks) +``` + +Each FAIL line names the case (``), the measured error and the +bound, e.g. `nfft_1d_50_50.txt … trafo_direct -> FAIL 5.7e+14 ( 1.07e-14)`. `perf-net.py` +guesses the **suite** from the filename prefix (`nfft` — `trafo_direct` is NFFT-specific, so +not `nfct`/`nfst`/`util`); that suite is the one to run in the inner loop. The detailed +machine-readable report is `tests/CUnitAutomated-Results.xml` (and `…_threads-Results.xml`). + +> **HARD GATE — no failing test ⇒ no coverage ⇒ stop.** If the injected fault leaves +> *every* test green (`perf-net.py check` reports GREEN, zero `-> FAIL`), the target is **not covered** by the +> suite. This is a blocking condition, not a green light: you cannot safely optimize +> code whose behaviour no test pins. Either add a test that fails on the fault first +> (then resume) — see [extending-tests.md](extending-tests.md) for how to add an online +> or refgen case — or stop and report the coverage gap. Never optimize an uncovered +> region. (Try a *more destructive* fault before concluding there is no coverage — a +> too-subtle edit may stay within tolerance.) + +**A4. Revert and re-confirm green.** Restore the exact original, rebuild, re-run; +`perf-net.py check` must report **GREEN** (exit 0) and `git diff` must be empty: + +```bash +cmake --build build-cmake -j >/dev/null +build-cmake/tests/checkall 2>&1 | uv run python $SCR/perf-net.py check - # GREEN, exit 0 +git diff --quiet && echo "diff empty" +``` + +You now know the precise tests that guard this region. + +**Seed the risk table.** While the net is fresh, note what it visibly does *not* cover +near the target — sizes, dims, node/coefficient distributions. Those gaps are the first +entries the [risk assessment](risk-assessment.md) carries forward into Phase D (the +[rounding-error analysis](phase-d-error-analysis.md)) and the inner loop. + +## Deliverables (exit criteria) + +Fill [`../templates/phase-b-correctness-net.md`](../templates/phase-b-correctness-net.md) +in `.perfeng/`. It has two outcomes — record exactly one: + +- **Net pinned ✅** — the fault flipped ≥1 case. The doc records: the injected fault + (saved verbatim as `artifacts/fault.diff`); the resulting net as a **Correctness + net** table (canonical format — see [deliverables.md](deliverables.md#canonical-formats)) + *with the suite to run in the inner loop and the net size*; and the revert + confirmation (`git diff` empty, suite green again). +- **Blocked ⛔** — no test failed even under a more destructive fault ⇒ region + uncovered. The doc is a **blocked report**: it documents the coverage gap and states + that the loop STOPS here — no Phase C/D/E. `artifacts/fault.diff` still records the + fault tried. + +*Deliverable = exit gate:* Phase B is not exitable until `phase-b-correctness-net.md` +and `artifacts/fault.diff` exist AND the tracker Phase B row is flipped — `✅` (net +pinned, proceed to C) or `⛔` (blocked, stop). On `⛔` the run ends here: set the +tracker header **Status** = `reverted` and write the human report +`summary.html` (`` = `fail`, documenting the coverage gap) — a blocked run +still gets a reviewer-facing report, then [conclude](phase-g-conclude.md) (there is no +optimization to land, so a PR is optional — offer one only if a permanent test closed the gap). diff --git a/.claude/skills/nfft-perf-eng/details/phase-c-performance-metric.md b/.claude/skills/nfft-perf-eng/details/phase-c-performance-metric.md new file mode 100644 index 00000000..f9b1fbdd --- /dev/null +++ b/.claude/skills/nfft-perf-eng/details/phase-c-performance-metric.md @@ -0,0 +1,73 @@ +# Phase C — pin the performance metric + +*[← Overview & map](../REFERENCE.md) · Prev: [Phase B — correctness net](phase-b-correctness-net.md) · Next: [Phase D — rounding-error analysis](phase-d-error-analysis.md)* + +**B1. Measure a baseline for the target's case(s)** with the walltime binary from +Phase A (run directly; it writes a per-case JSON): + +```bash +CODSPEED_PROFILE_FOLDER=/tmp/b build-cmake/benchmarks/bench_nfft_direct \ + --benchmark_filter='nfft_forward_direct_1d.*' +# /tmp/b/results/.json → per case: stats.median_ns (use the median, it is robust) +``` + +**B2. Inject a slowdown** in the target — e.g. wrap its body in a `for` loop that +repeats the work N times. Keep results correct (so you isolate cost, not correctness). + +> **Make the repeated work un-foldable — the build is `-O3 -ffast-math`.** A naive ×N +> wrapper *inside* the innermost body (recomputing the same value with the same arguments, +> then discarding all but the last) is silently removed by the optimizer: pure functions +> like `COS`/`SIN` with identical arguments are **common-subexpression-eliminated** to one +> call, and an overwritten accumulator is **dead-store-eliminated** — so the benchmark +> barely moves and the gate looks (wrongly) un-pinnable. Force work the optimizer cannot +> fold: wrap a **whole loop pass** (two full passes over the k-loop, which the compiler does +> not fuse), or make every repetition *observable* (e.g. accumulate all N copies into the +> live result and divide by N — `acc += term; … result = acc/N`). Verify the slowdown +> actually landed (the metric ≈doubles for a ×2) before trusting a "no movement" reading. +> A consistent few-percent rise across *all* sizes is the tell-tale of a partly-folded +> slowdown — not a pinned metric. + +**B3. Re-run and compare medians against the baseline.** Whichever benchmark cases' +`median_ns` rise clearly (well beyond `stdev_ns`) are the ones that exercise the +target — your performance metric. A case whose median is unchanged does **not** touch +the target. Collate both runs (`jq -s '[.[].benchmarks[]]' …/results/*.json > a.json`) +and let [`scripts/perf-bench.py compare`](../scripts/perf-bench.py) do the comparison +deterministically — under the *deliberate* slowdown the cases it flags `❌ regressed` +(`--base` = pre-slowdown, `--final` = slowed) are exactly your metric; here a "regression" +is the goal, so a non-zero exit confirms the metric moved. + +> **HARD GATE — no obtainable metric ⇒ no progress.** The agent must be able to +> produce a *concrete, comparable* **walltime** number (`median_ns`) for the target. +> Two ways this gate fails: +> (a) the benchmark tooling cannot run at all (see [Tooling status](tooling-status.md)), +> or (b) the tooling runs but the deliberate slowdown moves *no* benchmark — meaning +> no benchmark exercises the target, so a real improvement would be equally +> invisible. Either way, **stop**: add a benchmark that covers the target (see +> [caveats](caveats.md)) or pick a different target. Optimization without a metric is +> unverifiable and must not proceed. + +**B4. Revert** the slowdown; `git diff` must be empty. + +## Deliverables (exit criteria) + +Fill [`../templates/phase-c-performance-metric.md`](../templates/phase-c-performance-metric.md) +in `.perfeng/` and save +the injected slowdown verbatim as `artifacts/slowdown.diff`. Two outcomes — record +exactly one: + +- **Metric pinned ✅** — the doc records: the target's baseline measurement as a + [**Benchmark snapshot**](deliverables.md#canonical-formats) table (canonical format + — do not re-define columns); the injected slowdown (saved as `artifacts/slowdown.diff`); + the metric case(s) whose `median_ns` rose clearly, with before/after numbers; and + the revert confirmation (`git diff` empty). Flip the tracker Phase C row to ✅. +- **Blocked ⛔** — no benchmark moves *or* tooling can't measure (the HARD GATE, + branches (a)/(b) above). The deliverable is a **blocked report**: the doc explains + the no-metric verdict (uncovered target / tooling can't run), `artifacts/slowdown.diff` + is the slowdown that moved nothing. Set the tracker Phase C row to ⛔ — **the loop + STOPS**, do not proceed (no Phase D/E/F). + +*Deliverable = exit gate:* Phase C is not exitable until `phase-c-performance-metric.md` +exists, `artifacts/slowdown.diff` is saved, and the tracker row is flipped — ✅ (metric +pinned) ⇒ continue to D, or ⛔ (blocked) ⇒ the run ends here: like a Phase B block, set +the tracker header **Status** = `reverted`, write the human report +`summary.html` (`` = `fail`), then [conclude](phase-g-conclude.md) and report. diff --git a/.claude/skills/nfft-perf-eng/details/phase-d-error-analysis.md b/.claude/skills/nfft-perf-eng/details/phase-d-error-analysis.md new file mode 100644 index 00000000..02d08ceb --- /dev/null +++ b/.claude/skills/nfft-perf-eng/details/phase-d-error-analysis.md @@ -0,0 +1,174 @@ +# Phase D — rounding-error analysis (set the accuracy objective) + +*[← Overview & map](../REFERENCE.md) · Prev: [Phase C — performance metric](phase-c-performance-metric.md) · Next: [Phase E — inner loop](phase-e-inner-loop.md)* + +The two gates are green: the target is **covered** (B) and **measurable** (C). Before +changing a single line **for speed**, scrutinise the target's *numerical behaviour* with a +meticulous **rounding-error analysis under the standard model**. The optimization is then +held to a sharper standard than "tests still pass": a faster transform that is also *less +accurate* is usually not progress, and — read the other way — the analysis often reveals an +**avoidable** error source whose removal makes the code both **more accurate and faster**. + +This phase makes no optimization edits (the tree is clean again after C's revert). It only +*reads, derives, and measures*, and emits a verdict that **sets the objective** of the +inner loop: + +- **`clean`** — the target's error is at (or provably near) the best achievable floor for + what it computes. Objective for Phase E: **pure performance**, with *no accuracy + regression* — and now the derived bound, not just the pass/fail tests, is the yardstick. +- **`improve-first`** — the analysis and/or the baseline test errors expose **error beyond + the floor**: an *avoidable* dependence on the transform size `N`/`M`, catastrophic + cancellation, a poorly conditioned intermediate, an unstable summation/recurrence. + Objective for Phase E: **first root out the avoidable error, then tune speed.** Often the + same restructuring does both ("often enough you win on both metrics"). + +## D1. The standard model + +State and use the standard model of floating-point arithmetic (unit roundoff `u`): + +> for `∘ ∈ {+,−,×,÷}`, `fl(a ∘ b) = (a ∘ b)(1 + δ)`, `|δ| ≤ u`, and likewise +> `fl(√a) = √a (1 + δ)`. For a length-`n` accumulation the errors compose into the +> Higham constant `γₙ = n·u / (1 − n·u) ≈ n·u`. + +`u` per build precision (the library's `MANT_DIG` branches): float `2⁻²⁴`, double `2⁻⁵³`, +long double `2⁻⁶⁴` (80-bit) / `2⁻¹¹³` (quad). The analysis is **precision-agnostic** in +form (a `γ`-bound in `n` and `u`) and instantiates per precision — the same reason the loop +measures all three (see [precision-matrix](precision-matrix.md)). + +## D2. Derive the target's error bound (theory) + +Write the target as a sequence of floating-point operations and bound the forward error, +**keeping the dependence on `N` and `M` explicit** — that dependence is the whole point. + +Worked example — the direct NDFT `f_j = Σ_{k} f̂_k e^{-2πi k x_j}` (`X(trafo_direct)`): + +1. Each summand is a product `f̂_k · w` with `w = e^{-2πi k x_j}` formed by `COS`/`SIN` + (or, after an optimization, a recurrence). Bound the **per-summand** error first: the + transcendental's relative error, plus argument-reduction error in `k·x_j`, plus the + complex multiply. +2. The length-`N` **sequential (recursive) summation** contributes the Higham factor + `γ_{N-1}`. Give **both** bounds, because they differ sharply for this transform's data: + - **Worst case (guarantee):** `E∞/‖f̂‖₁ ≤ γ_{N-1} + c·u ≈ (N + c)·u` (measured as the suite + does — see [`test-methodology.md`](../../../../docs/agents/test-methodology.md); the + normalisation works because `|e^{-2πi k x_j}| = 1`, so `Σ|f̂_k·w_k| = ‖f̂‖₁`). + - **Expected (random signs):** the test inputs have random signs, so the rounding errors + partially cancel — the RMS error grows like `~√N·u`, the full `N·u` being a pessimistic + envelope rarely realised. `√N·u` is the *dominant* term random data exhibits. (In a + measured trend the fitted exponent often lands a little **above** 0.5 — ≈0.6–0.7 — because + a sub-dominant working-precision term, e.g. forming the phase `2π(k−N/2)x_j` to absolute + error `~N·u`, lifts it; this shows up in the *baseline* too and is not an alarm. See D3.) +3. Compare with the suite's **bound**: the direct transforms use a *flat* round-off floor + `bound = C·ε` (`C` = 48 NFFT / 120 NFCT / 130 NFST). Under the **√N** regime the test data + exhibits, a flat constant `C` is robust across a wide size range — *not* a fragile + coincidence. Only the **worst-case** `N·u` envelope would eventually exceed `C·ε`; flag + that as a *worst-case-only* `size-dependent` note (D3 decides from the measured trend + whether it is a real risk) — never raise the alarm on the `N·u` envelope alone when the + data sits on `√N` ([risk-assessment](risk-assessment.md)). + +Separate **inherent** from **avoidable**: summing `N` terms is fundamentally `O(N·u)` +(mitigable by compensated summation, never zero); an implementation artifact — a recurrence +that compounds to *worse* than `N·u`, a cancelling reformulation, an un-reduced large +argument — is **avoidable** and is exactly what `improve-first` targets. + +## D2′. Calibrate to where the error matters (don't minimise blindly) + +The objective is **not** "drive the target's error to the smallest possible". The right +amount of accuracy depends on the target's role in the larger computation — situate it: + +- **Downstream-limited.** If a *later* stage deliberately approximates (a truncated window, + a fixed cutoff `m`), the end-to-end accuracy is capped there. Pushing the target below that + cap with heavy machinery (compensated, or **doubly**-compensated summation, exact products) + buys nothing the user sees — and may cost speed. Match the target's error to the cap, no + finer. +- **Leveraged / amortised.** If the target's result is *reused* — above all a **precomputation** + feeding many transforms — its error propagates into every consumer, and the cost of computing + it more accurately is paid **once**. Here a smaller error is highly desirable even at a real + one-off cost: weigh per-call vs one-off explicitly. +- **Terminal.** If the target *is* (part of) the returned result, its error is the user's + error; the derived bound is the thing to protect. + +State which case the target is in — it sets *how hard* `improve-first` should push, and it is +what makes an accuracy/speed trade-off a **reasoned** call rather than a reflex. + +## D3. Cross-check against the baseline (measurement) + +Theory predicts an *order of growth*; the Phase-A baseline already measured the *actual* +error. Reconcile them: + +- The Phase-A `ctest` logs (`artifacts/baseline-tests-{d,f,l}.log`) print each case's + measured error and its bound. Tabulate **measured error vs `N` (and `M`)** for the + target's cases, per precision. Expect the trend dominated by **`√N`** (random-sign + cancellation), well below the worst-case `N·u` envelope — that is the *healthy* baseline. + A fitted exponent a little above 0.5 (≈0.6–0.7) is normal — the sub-dominant + working-precision phase term (above), present in the baseline too; the `trafo-direct` + reference run measured + ≈0.62–0.68 in *both* arms. What would be alarming is the **optimized** trend growing + *faster than the baseline* (an avoidable source the theory pinned, or the worst case being + realised) — i.e. a non-trivial Δp, not the absolute `p`. +- When the regime is in question, run the **differential trend study in all three precisions** + (float and long double can diverge — float cancels/overflows sooner, so a double-only trend + hides a precision-specific change): a geometric sweep of `N` against an oracle **strictly more + precise than the build** — `tests.refgen.transforms` at high `mpmath` `dps`, *not* a + same-precision C reference (a quad C oracle cannot validate the long-double build; see + [extending-tests.md](extending-tests.md#differential-trend-analysis-the-strongest-settle-for-accuracy-for-speed)). + Fit `error ≈ C·Nᵖ` with [`scripts/perf-trend.py`](../scripts/perf-trend.py) and **compare the + exponents** per precision — the verdict is `Δp = p_opt − p_base` (≈0 ⇒ preserved), not the + absolute `p`. The data + lands as `artifacts/trend-{baseline,optimized}-{d,f,l}.dat` and is charted tabbed (one per + precision) in `summary.html` by [`scripts/perf-summary.py`](../scripts/perf-summary.py). This is the empirical half of the analysis and the same + tool Phase E/the risk loop use to settle accuracy-for-speed risks — see + [extending-tests.md](extending-tests.md#differential-trend-analysis-the-strongest-settle-for-accuracy-for-speed). + +## D4. Verdict — classify and set the objective + +Conclude with one of `clean` / `improve-first`, the derived bound (per precision), the +measured-vs-derived reconciliation, the **role classification** from D2′ (downstream-limited +/ leveraged / terminal), and — for `improve-first` — the **named avoidable source** and how +the optimization should remove it. This verdict is what [Phase E](phase-e-inner-loop.md) +optimizes toward and what [Phase F](phase-f-exit-gate.md) checks was honoured. It also +**seeds the risk table** with the size/precision-dependence the analysis exposed, so the +inner loop carries it forward ([risk-assessment](risk-assessment.md)). + +> The honest outcome may be **`clean`: the current error behaviour is the best achievable** +> for what the target computes — that is a valid, common result, not a failure to find work. +> Say so plainly and optimize for speed alone (without regressing accuracy). + +**When accuracy and speed genuinely conflict, escalate — don't decide silently.** The happy +path is one restructuring that improves both. But a *real* trade-off exists — removing the +avoidable error **costs** speed with no joint win (e.g. doubly-compensated summation on a +target whose accuracy is anyway downstream-limited, or a more accurate but slower +precomputation amortised differently than you'd guess). That is a **human-judgement call**, +not the agent's to make unilaterally: the loop is autonomous, but this decision is the +reviewer's. Pick a documented default (the role classification from D2′ is the argument — +protect leveraged/terminal accuracy, don't over-invest in downstream-limited), record the +trade-off and the alternative explicitly in `error-analysis.html`, and carry it to +[Phase F](phase-f-exit-gate.md) so the run surfaces it (a `partial` outcome with the trade-off +laid out) rather than burying a chosen point on the accuracy/speed curve. + +## Deliverables (exit criteria) + +Two deliverables, both in `.perfeng/`: + +- **`error-analysis.html`** — the **primary** deliverable: a self-contained, MathJax-rendered + document carrying the full standard-model derivation (the math the reviewer must follow), + the measured-vs-derived reconciliation, and the verdict. Copy + [`../templates/error-analysis.html`](../templates/error-analysis.html) and fill it; it is + self-contained (MathJax via CDN ` + + + + + +

+ +
Rounding-error analysis · NFFT3 · Phase D
+

X(target)() — error under the standard model

+

+ +
+ <✅ clean — error at the best achievable floor | ⚠ improve-first — avoidable error source> — + +
+ TargetX(target)() + Metric\(E_\infty/\lVert \hat f\rVert_1\) — see phase-d-error-analysis.md + Baseline commit + TrackerREADME.md +
+
+ +

1 · The standard model

+

All bounds use the standard model of floating-point arithmetic with unit roundoff + \(u\) (\(u=2^{-24}\) float, \(2^{-53}\) double, \(2^{-64}\) / \(2^{-113}\) long + double — the library's MANT_DIG branches):

+ \[ \mathrm{fl}(a \circ b) = (a \circ b)(1+\delta), \qquad |\delta| \le u, \qquad \circ \in \{+,-,\times,\div\}. \] +

For a length-\(n\) accumulation the per-operation errors compose into the Higham + constant

+ \[ \gamma_n = \frac{n\,u}{1 - n\,u}, \qquad \gamma_n \le 1.01\,n\,u \ \ (n u \le 0.01), \qquad \gamma_n \approx n\,u . \] + +

2 · The computation as floating-point operations

+

+
<verbatim target code or the hot expression — kernel/.../file.c:NN>
+ + \[ f_j = \sum_{k} \hat f_k \, e^{-2\pi i k x_j}, \qquad + \hat f_j := \mathrm{fl}\!\left(\sum_{k} \hat f_k\, \mathrm{fl}\!\big(e^{-2\pi i k x_j}\big)\right). \] +

+ +

3 · Forward-error bound (keep \(N\), \(M\) explicit)

+

+ + \[ \underbrace{\frac{|\hat f_j - f_j|}{\lVert \hat f\rVert_1} \;\le\; \gamma_{N-1} + c\,u \;\approx\; (N + c)\,u}_{\text{worst case (guarantee)}}, + \qquad + \underbrace{\mathbb{E}\!\left[\frac{|\hat f_j - f_j|}{\lVert \hat f\rVert_1}\right] \;\sim\; \sqrt{N}\,u}_{\text{expected, random signs}} . \] +

where \(c\) collects the per-summand transcendental and multiply terms (the + normalisation works because \(|e^{-2\pi i k x_j}| = 1\), so \(\sum|\hat f_k w_k| = + \lVert\hat f\rVert_1\)). The worst case is the deterministic envelope; with random-sign + coefficients the errors partly cancel and the expected error grows like \(\sqrt N\,u\). +

+ + +

4 · Measured vs derived (baseline reconciliation)

+

artifacts/baseline-tests-{d,f,l}.log. Expect + \(p\approx 0.5\) (\(\sqrt N\)) for random-sign data — \(p\to 1\) signals the worst case + is realised or an avoidable source. If a differential trend study was run, cite the + fitted exponent \(p\) in \(\text{error}\approx C N^{p}\) from + scripts/perf-trend.py.>

+ + + + + +
prec\(N\) (\(M\))measured errorboundtrend \(p\) (expect \(\approx 0.5\))
d<…><…><…><\(p\approx 0.5\) — \(\sqrt N\)>
f<…><…><…><…>
l<…><…><…><…>
+

+ +

5 · Verdict & objective for the inner loop

+

Target role:

+

clean or improve-first and why, per precision. For improve-first, + NAME the avoidable source and how the optimization removes it — and note where both + accuracy and speed can improve together. For clean, say plainly that the current + behaviour is the best achievable for what the target computes, so Phase E optimizes for + speed alone without regressing the derived bound.>

+ + +

6 · Risks this analysis surfaces

+ + + + + +
RiskCategoryStateEvidence / disposition
+ + +
+ + diff --git a/.claude/skills/nfft-perf-eng/templates/phase-a-baseline.md b/.claude/skills/nfft-perf-eng/templates/phase-a-baseline.md new file mode 100644 index 00000000..854a7852 --- /dev/null +++ b/.claude/skills/nfft-perf-eng/templates/phase-a-baseline.md @@ -0,0 +1,22 @@ +# Phase A — baseline (the exit reference) + +- **Build:** `` (three trees: `build-cmake{,-f,-l}`) +- **Baseline commit:** +- **Tests (per precision):** double ✅ · float ✅ · long double ✅ — `artifacts/baseline-tests-{d,f,l}.log` +- **Benchmark cases captured:** × 3 precisions — `artifacts/baseline-bench-{d,f,l}.json` + +## Benchmark snapshot + + +| prec | case | median_ns | stdev_ns | rounds | +|------|----------------------------|-----------|----------|--------| +| d | | <123456> | <789> | <50> | +| f | | <98765> | <654> | <50> | +| l | | <210987> | <912> | <50> | + +## Notes + + +— diff --git a/.claude/skills/nfft-perf-eng/templates/phase-b-correctness-net.md b/.claude/skills/nfft-perf-eng/templates/phase-b-correctness-net.md new file mode 100644 index 00000000..899fde47 --- /dev/null +++ b/.claude/skills/nfft-perf-eng/templates/phase-b-correctness-net.md @@ -0,0 +1,35 @@ +# Phase B — correctness net + +- **Outcome:** +- **Target:** `` — +- **Fault injected:** — `artifacts/fault.diff` +- **Suite to run in inner loop:** `` +- **Net size:** cases + +## Correctness net + + +| suite | case | error | bound | +|--------|------------------------------------|-----------|------------| +| | | <5.7e+14> | <1.07e-14> | + +## Revert confirmation + +- suite green again: · `git diff` empty: + + + + diff --git a/.claude/skills/nfft-perf-eng/templates/phase-c-performance-metric.md b/.claude/skills/nfft-perf-eng/templates/phase-c-performance-metric.md new file mode 100644 index 00000000..05f55b55 --- /dev/null +++ b/.claude/skills/nfft-perf-eng/templates/phase-c-performance-metric.md @@ -0,0 +1,40 @@ +# Phase C — performance metric + +- **Outcome:** +- **Benchmark filter:** `<--benchmark_filter=…>` +- **Slowdown injected:** — `artifacts/slowdown.diff` + +## Target baseline + + +| case | median_ns | stdev_ns | rounds | +|----------|-----------|----------|--------| +| | | | | + +## Metric — cases that moved under the slowdown + + +| case | median_ns base | median_ns slowed | moves the target? | +|--------|----------------|------------------|-------------------| +| | | | | + +## Revert confirmation + +- `git diff` empty: + + + + diff --git a/.claude/skills/nfft-perf-eng/templates/phase-d-error-analysis.md b/.claude/skills/nfft-perf-eng/templates/phase-d-error-analysis.md new file mode 100644 index 00000000..01865041 --- /dev/null +++ b/.claude/skills/nfft-perf-eng/templates/phase-d-error-analysis.md @@ -0,0 +1,52 @@ +# Phase D — rounding-error analysis (verdict & objective) + + + +- **Verdict:** +- **Target role:** +- **Full analysis:** [error-analysis.html](error-analysis.html) +- **Standard model `u`:** float `2⁻²⁴` · double `2⁻⁵³` · long double `2⁻⁶⁴`/`2⁻¹¹³` + +## Derived bound + + +- **Worst case (guarantee):** +- **Expected (random signs):** <`~√N·u` — what random data exhibits; the N·u envelope is loose> +- **Suite bound vs derived:** ) robust under the √N regime; only the worst-case `N·u` envelope breaches, at N≈<…> — worst-case-only> + +## Measured vs derived (baseline reconciliation) + + +| prec | N (M) | measured error | bound | trend p (expect ≈0.5) | +|------|-------|----------------|-------|-----------------------| +| d | <…> | <…> | <…> | | + +## Objective for Phase E + + +<...> + + +- **Trade-off:** costs % speed — default , alternative ; reviewer to weigh> + +## Seeded risks (roll into the Phase-E/F risk table) + + +| risk | category | state | evidence / disposition | +|----------------------------------------|----------------|----------|---------------------------------------| +| | | | | + + diff --git a/.claude/skills/nfft-perf-eng/templates/phase-e-inner-loop.md b/.claude/skills/nfft-perf-eng/templates/phase-e-inner-loop.md new file mode 100644 index 00000000..81700a68 --- /dev/null +++ b/.claude/skills/nfft-perf-eng/templates/phase-e-inner-loop.md @@ -0,0 +1,38 @@ +# Phase E — inner loop (iteration journal) + +- **B-net:** `` ( cases) · **C-metric:** +- **Phase-D objective:** first, then tune speed> +- **Current change:** `artifacts/change.diff` + +## Iterations + + +| iter | change | net | metric median (ns) before→after | +|------|--------------------------------|--------------------------------------|----------------------------------| +| 1 | | | | +| 2 | | | <— (reverted)> | + +## Risks raised (rolls up to Phase F) + + +| risk | category | state | evidence / disposition | +|------------------------------------------|----------------|----------|------------------------------------------| +| | | | | + +## Exit state + +- B-net green at latest kept state, **all precisions**: +- metric median dropped beyond noise: +- Phase-D accuracy objective met: +- `artifacts/change.diff` matches latest kept state: +- risk table current; every material risk `proven`/`retired`/`accepted`/`residual`: + + diff --git a/.claude/skills/nfft-perf-eng/templates/phase-f-exit-gate.md b/.claude/skills/nfft-perf-eng/templates/phase-f-exit-gate.md new file mode 100644 index 00000000..9ad1131e --- /dev/null +++ b/.claude/skills/nfft-perf-eng/templates/phase-f-exit-gate.md @@ -0,0 +1,62 @@ +# Phase F — exit gate + +- **Verdict:** +- **Target metric:** : ns (<Δ%>) +- **Phase-D objective:** > — honoured: +- **Landed change:** `artifacts/change.diff` + +## Comparison (baseline → final, all cases × all precisions) + + +| prec | case | base median_ns | final median_ns | Δ% | threshold | verdict | +|------|--------------------------|----------------|-----------------|------|-----------|-----------| +| d | | <123456> | <95012> | <−23%> | 2%/3σ | <✅ faster> | +| f | | <98765> | <80120> | <−19%> | 2%/3σ | <✅ faster> | +| l | | <210987> | <165430> | <−22%> | 2%/3σ | <✅ faster> | + +## Exit conditions (all must hold, in float · double · long double) + +1. Full suite passes as in Phase A (incl. `checkall_threads`), every precision: — `artifacts/final-tests-{d,f,l}.log` +2. No benchmark regresses beyond noise (every case, every precision) — `perf-bench.py compare --taskdir .`: +3. Surviving regressions on untouched cases attributed (callgrind `I refs`: layout artifact vs real coupling): +4. `git diff` is only the intended change: +5. Phase-D accuracy objective honoured — `improve-first`: avoidable source removed, error no worse (trend study if order-of-growth); `clean`: derived bound not regressed: +6. Risk assessment complete & honest — every material risk (incl. Phase-D-seeded rows) `proven`/`retired`/`accepted`/`residual`; no `proven` material drop landed (fixed or reverted); any unsettled `residual` → outcome `partial`: + +## Verdict + + +<...> + +## Risk assessment (consolidated risk table) + + +| risk | category | state | evidence / disposition | +|-----------------------------------------|----------------|----------|-----------------------------------------| +| | | | | +| | | | | + +Note: `proven` material accuracy drop is a hard no — it never reaches this consolidated table as a *landed* state (it was fixed or the run reverted). An unsettled `residual` material drop ⇒ this run's outcome is `partial`. + +## Close-out + + +- temporary probing tests removed (their findings kept in the risk table); permanent test additions noted as expected in the `git diff` +- tracker header **Status** = +- tracker **Outcome** one-liner filled +- tracker **Phase F** row flipped (`✅` complete / `⛔` reverted) + Exit signal set +- **charts generated** — `perf-summary.py charts --taskdir .` → `artifacts/chart-speedup-{d,f,l}.svg` (+ `chart-trend-{d,f,l}.svg` if a trend study ran) +- **`summary.html` written** (from `../templates/summary.html`) — the human report; + presents every phase's result + numbers, embeds the required charts, links every deliverable + artifact; + set `` to `ok` (completed clean) · `partial` (landed with caveats — an unsettled `residual` material risk) · `fail` (reverted / hard-gate block) +- **completeness verified** — `perf-summary.py check --taskdir .` exits 0 (nothing orphaned; required charts linked) +- **concluded (Phase G)** — `perf-conclude.sh squash -m "…"` squashed the run to one commit; then (opt-in) pushed, opened a `perf-eng`-labelled PR to `develop`, and `perf-conclude.sh package ` renamed `.perfeng/` → `.perfeng-pr-/` + zipped it, attached to the PR diff --git a/.claude/skills/nfft-perf-eng/templates/summary.html b/.claude/skills/nfft-perf-eng/templates/summary.html new file mode 100644 index 00000000..0aa4bdcc --- /dev/null +++ b/.claude/skills/nfft-perf-eng/templates/summary.html @@ -0,0 +1,316 @@ + + + + + + +Perf · <target> — optimization summary + + + + +
+ +
Performance optimization · NFFT3
+

X(target)()

+

+ + +
+ <✅ Completed | ⚠ Landed with caveats | ⛔ Reverted / blocked> — + +
+ TargetX(target)() + Metricwalltime (median_ns) — local + CI (CodSpeed Macro Runners) + Precisions + Baseline commit + Date + TrackerREADME.md · phase deliverables linked in the footer +
+
+ +

1 · The target — what we optimized and why

+

+ +
<verbatim target code or the hot snippet — kernel/.../file.c:NN>
+

+ +

2 · The idea — what to change, and the reasoning

+

+ + +

3 · The process — idea to verdict

+ + + + + + + + + +
PhaseWhat happened
A · baseline
B · correctness netnfft).>
C · performance metric
D · error analysisclean | improve-first: <source>; the accuracy objective it set. See error-analysis.html.>
E · inner loop
F · exit gate
+ +

4 · Performance results

+ + +
+ + + +
+
double — % faster per case (base→final)
+
float — % faster per case (base→final)
+
long double — % faster per case (base→final)
+
+
+

% faster per case (green faster · red regressed · grey within noise, threshold + max(3σ, 2%)). One tab per precision; generated by perf-summary.py charts from + the baseline/final benchmark JSON (linked below).

+ + + + + + + + +
preccasebefore (ns)after (ns)Δnote
d<123456><95012><−23%>
f<98765><80120><−19%>
l<210987><165430><−22%>
+

+ +
    +
  • <✅/❌> Full suite passes as Phase A (d/f/l, incl. checkall_threads)
  • +
  • <✅/❌> No benchmark regresses beyond noise (every case, every precision)
  • +
  • <✅/N/A> Surviving regressions on untouched cases attributed (callgrind I refs: layout vs real coupling)
  • +
  • <✅/❌> git diff is only the intended change
  • +
  • <✅/❌> Phase-D accuracy objective honoured (avoidable source removed / derived bound held)
  • +
  • <✅/❌> Risk assessment complete & honest — no proven material drop landed
  • +
+ +

5 · Accuracy & error behaviour

+ +

Verdict: · role: . + Full derivation: error-analysis.html.

+ +
+ + + +
+
double — error vs N (log-log), baseline vs optimized
+
float — error vs N (log-log), baseline vs optimized
+
long double — error vs N (log-log), baseline vs optimized
+
+
+

Error growth vs N (log-log), baseline vs optimized, fitted + exponents p (expect ≈0.5, √N); one tab per precision. Generated by + perf-summary.py charts from trend-{baseline,optimized}-{d,f,l}.dat + (linked below).

+ + + +
preccase (N, M)error beforeerror afterboundtrend p
d<…><…><…><…><≈0.5>
+ +

6 · Risk assessment — what this result might still hide

+ +

+ + + + + + +
RiskCategoryStateEvidence / disposition
+

State: proven = a material accuracy drop, a hard no (fixed or + reverted — never lands) · retired = a check that would expose it now passes · + accepted = shown within the documented error bound, or a non-accuracy concern + reasoned through · residual = a suspected material drop no constructible check + could settle. A residual material risk doesn't block landing but makes this run + partial, not ok — "we didn't look" and "we looked and couldn't settle it" are + different, and only the second may land (for the reviewer to weigh).

+ + +

7 · Why it stalled (failure analysis)

+

+
    +
  • Known:
  • +
  • Hypothesis:
  • +
  • To unblock:
  • +
+ + +

8 · Patterns

+
+
+

Used

+
+
+
+

Discovered

+
+
+
+ +
+

9 · Carry-forward notes

+

Reusable when optimizing other targets in this codebase — accumulate the + knowledge so the next run starts ahead.

+ +
+ + + +
+ + diff --git a/.claude/skills/nfft-perf-eng/templates/tracker.md b/.claude/skills/nfft-perf-eng/templates/tracker.md new file mode 100644 index 00000000..473dfa67 --- /dev/null +++ b/.claude/skills/nfft-perf-eng/templates/tracker.md @@ -0,0 +1,30 @@ +# Perf: + +- **Target:** `` — +- **Baseline commit:** +- **Status:** in-progress +- **Outcome:** — + +## Phase status + +| Phase | Status | Deliverable | Exit signal | +|-------|--------|-------------|-------------| +| Phase A — baseline | 🔄 | [phase-a-baseline.md](phase-a-baseline.md) | — | +| Phase B — net | ⬜ | — | — | +| Phase C — metric | ⬜ | — | — | +| Phase D — error analysis | ⬜ | — | — | +| Phase E — inner | ⬜ | — | — | +| Phase F — exit | ⬜ | — | — | +| Phase G — conclude | ⬜ | — (squash + PR) | — | + + + +## Log + +- — created; baseline commit . diff --git a/.gitignore b/.gitignore index b3980494..ea46b250 100644 --- a/.gitignore +++ b/.gitignore @@ -37,6 +37,14 @@ # Log files *.log +# nfft-perf-eng deliverables — never committed. The single static working directory the +# skill writes every run into; its record ships as the run's squashed commit + PR + the +# zip attached to that PR, not in the repo tree. (.perfeng.bak = perf-init's one-level undo; +# .perfeng-pr- = the run renamed at conclude, and where a downloaded results zip unpacks.) +/.perfeng/ +/.perfeng.bak/ +/.perfeng-pr-*/ + # Autosave files *.asv @@ -159,7 +167,7 @@ matlab/tests/nfsftUnitTests.m matlab/tests/check_*_matlab.output windows-build-dll*/* linux-build-mex/* -nfft-* +/nfft-* .vscode/* .venv/ @@ -172,4 +180,4 @@ build-cmake*/ .DS_Store # JetBrains IDEs -.idea +.idea \ No newline at end of file diff --git a/kernel/nfct/nfct.c b/kernel/nfct/nfct.c index cdbe1070..d02d3bf3 100644 --- a/kernel/nfct/nfct.c +++ b/kernel/nfct/nfct.c @@ -64,6 +64,21 @@ static inline INT intprod(const INT *vec, const INT a, const INT d) #define NODE(p,r) (ths->x[(p) * ths->d + (r)]) +/* Block size for the phase recurrence in the direct transforms */ +#define NFFT_DIRECT_RECURRENCE_BLOCK 32 + +/* The (co)sine value is a part of the complex phase exp(+i 2pi (k+OFFSET) x): NDCT reads the + * real part (cos), NDST reads the imaginary part (sin). */ +#define BASEPART(w) CREAL(w) + +/* Accurate phase for exp(+i 2pi k x): reduce k*x modulo 1 into ~[-1/2,1/2) so COS/SIN see a + * small argument, error does not grow with N. Requires FMA single-rounding semantics. */ +static inline R X(reduced_omega)(const R k, const R x) +{ + const R n = RINT(k * x); // Nearest integer to k * x. + return K2PI * FFMA(k, x, -n); // Calculate k * x - n with a single rounding, then multiply 2 * pi. +} + #define MACRO_with_FG_PSI fg_psi[t][lj[t]] #define MACRO_with_PRE_PSI ths->psi[(j * ths->d + t) * (2 * ths->m + 2) + lj[t]] #define MACRO_without_PRE_PSI PHI((2 * NN(ths->n[t])), ((ths->x[(j) * ths->d + t]) \ @@ -94,18 +109,32 @@ void X(trafo_direct)(const X(plan) *ths) if (ths->d == 1) { /* specialize for univariate case, rationale: faster */ + const INT B = NFFT_DIRECT_RECURRENCE_BLOCK; INT j; #ifdef _OPENMP #pragma omp parallel for default(shared) private(j) #endif for (j = 0; j < ths->M_total; j++) { - INT k_L; - for (k_L = 0; k_L < ths->N_total; k_L++) + R v = K(0.0); + const R x = ths->x[j]; + const R dphi = K2PI * x; /* |dphi| <= pi: accurate without reduction */ + const C dw = COS(dphi) + II * SIN(dphi); /* per-step phase factor exp(+i 2pi x) */ + INT k_L = 0; + while (k_L < ths->N_total) { - R omega = K2PI * ((R)(k_L + OFFSET)) * ths->x[j]; - f[j] += f_hat[k_L] * BASE(omega); + /* Accurate seed exp(+i 2pi (k_L + OFFSET) x), then recur within the block. */ + const R omega = X(reduced_omega)((R)(k_L + OFFSET), x); + C w = COS(omega) + II * SIN(omega); + INT kend = k_L + B; if (kend > ths->N_total) kend = ths->N_total; + for (; k_L < kend; k_L++) + { + v += f_hat[k_L] * BASEPART(w); + w *= dw; + } } + + f[j] = v; } } else @@ -156,7 +185,49 @@ void X(adjoint_direct)(const X(plan) *ths) if (ths->d == 1) { /* specialize for univariate case, rationale: faster */ + const INT B = NFFT_DIRECT_RECURRENCE_BLOCK; #ifdef _OPENMP + if (ths->N_total > B) + { + /* Give each thread a disjoint, contiguous range of frequencies [klo,khi) (so the + * f_hat[k] writes are race-free) and run the phase recurrence within it, re-seeded + * every B steps. Cap the team at N/B threads so every thread owns at least one full + * block (>= B frequencies): with sub-block ranges the per-block seed (a COS+SIN pair) + * is un-amortised and the recurrence loses to the plain per-k loop. The cap depends + * only on N and B, so it removes any dependence on the machine's thread count. */ + const int nt_cap = (int)(ths->N_total / B); + #pragma omp parallel default(shared) num_threads(nt_cap) + { + const int nt = omp_get_num_threads(); + const int tid = omp_get_thread_num(); + const INT klo = (INT)(((long long)ths->N_total * tid) / nt); + const INT khi = (INT)(((long long)ths->N_total * (tid + 1)) / nt); + INT j; + for (j = 0; j < ths->M_total; j++) + { + const R x = ths->x[j]; + const R dphi = K2PI * x; + const C dw = COS(dphi) + II * SIN(dphi); + INT k_L = klo; + while (k_L < khi) + { + const R omega = X(reduced_omega)((R)(k_L + OFFSET), x); + C w = COS(omega) + II * SIN(omega); + INT kend = k_L + B; if (kend > khi) kend = khi; + for (; k_L < kend; k_L++) + { + f_hat[k_L] += f[j] * BASEPART(w); + w *= dw; + } + } + } + } + } + else + { + /* N <= B: the recurrence spans at most one block per thread-range, so its per-block + * seed/setup costs more than it saves once threaded. Use the plain per-k + * parallelisation. At these tiny N the per-entry phase error is in check. */ INT k_L; #pragma omp parallel for default(shared) private(k_L) for (k_L = 0; k_L < ths->N_total; k_L++) @@ -168,15 +239,43 @@ void X(adjoint_direct)(const X(plan) *ths) f_hat[k_L] += f[j] * BASE(omega); } } + } #else INT j; - for (j = 0; j < ths->M_total; j++) + if (ths->N_total > B) { - INT k_L; - for (k_L = 0; k_L < ths->N_total; k_L++) + for (j = 0; j < ths->M_total; j++) { - R omega = K2PI * ((R)(k_L + OFFSET)) * ths->x[j]; - f_hat[k_L] += f[j] * BASE(omega); + const R x = ths->x[j]; + const R dphi = K2PI * x; + const C dw = COS(dphi) + II * SIN(dphi); + INT k_L = 0; + while (k_L < ths->N_total) + { + const R omega = X(reduced_omega)((R)(k_L + OFFSET), x); + C w = COS(omega) + II * SIN(omega); + INT kend = k_L + B; if (kend > ths->N_total) kend = ths->N_total; + for (; k_L < kend; k_L++) + { + f_hat[k_L] += f[j] * BASEPART(w); + w *= dw; + } + } + } + } + else + { + /* N <= B: the recurrence spans at most one block, so its per-j seed/setup outweighs + * the saved transcendentals (which the plain loop vectorises); the un-reduced + * argument's phase error is still in check at these tiny N. */ + for (j = 0; j < ths->M_total; j++) + { + INT k_L; + for (k_L = 0; k_L < ths->N_total; k_L++) + { + R omega = K2PI * ((R)(k_L + OFFSET)) * ths->x[j]; + f_hat[k_L] += f[j] * BASE(omega); + } } } #endif diff --git a/kernel/nfst/nfst.c b/kernel/nfst/nfst.c index df71f4a1..e0c70212 100644 --- a/kernel/nfst/nfst.c +++ b/kernel/nfst/nfst.c @@ -64,6 +64,21 @@ static inline INT intprod(const INT *vec, const INT a, const INT d) #define NODE(p,r) (ths->x[(p) * ths->d + (r)]) +/* Block size for the phase recurrence in the direct transforms */ +#define NFFT_DIRECT_RECURRENCE_BLOCK 32 + +/* The (co)sine value is a part of the complex phase exp(+i 2pi (k+OFFSET) x): NDCT reads the + * real part (cos), NDST reads the imaginary part (sin). */ +#define BASEPART(w) CIMAG(w) + +/* Accurate phase for exp(+i 2pi k x): reduce k*x modulo 1 into ~[-1/2,1/2) so COS/SIN see a + * small argument, error does not grow with N. Requires FMA single-rounding semantics. */ +static inline R X(reduced_omega)(const R k, const R x) +{ + const R n = RINT(k * x); // Nearest integer to k * x. + return K2PI * FFMA(k, x, -n); // Calculate k * x - n with a single rounding, then multiply 2 * pi. +} + #define MACRO_with_FG_PSI fg_psi[t][lj[t]] #define MACRO_with_PRE_PSI ths->psi[(j * ths->d + t) * (2 * ths->m + 2) + lj[t]] #define MACRO_without_PRE_PSI PHI((2 * NN(ths->n[t])), ((ths->x[(j) * ths->d + t]) \ @@ -94,18 +109,32 @@ void X(trafo_direct)(const X(plan) *ths) if (ths->d == 1) { /* specialize for univariate case, rationale: faster */ + const INT B = NFFT_DIRECT_RECURRENCE_BLOCK; INT j; #ifdef _OPENMP #pragma omp parallel for default(shared) private(j) #endif for (j = 0; j < ths->M_total; j++) { - INT k_L; - for (k_L = 0; k_L < ths->N_total; k_L++) + R v = K(0.0); + const R x = ths->x[j]; + const R dphi = K2PI * x; /* |dphi| <= pi: accurate without reduction */ + const C dw = COS(dphi) + II * SIN(dphi); /* per-step phase factor exp(+i 2pi x) */ + INT k_L = 0; + while (k_L < ths->N_total) { - R omega = K2PI * ((R)(k_L + OFFSET)) * ths->x[j]; - f[j] += f_hat[k_L] * BASE(omega); + /* Accurate seed exp(+i 2pi (k_L + OFFSET) x), then recur within the block. */ + const R omega = X(reduced_omega)((R)(k_L + OFFSET), x); + C w = COS(omega) + II * SIN(omega); + INT kend = k_L + B; if (kend > ths->N_total) kend = ths->N_total; + for (; k_L < kend; k_L++) + { + v += f_hat[k_L] * BASEPART(w); + w *= dw; + } } + + f[j] = v; } } else @@ -156,7 +185,49 @@ void X(adjoint_direct)(const X(plan) *ths) if (ths->d == 1) { /* specialize for univariate case, rationale: faster */ + const INT B = NFFT_DIRECT_RECURRENCE_BLOCK; #ifdef _OPENMP + if (ths->N_total > B) + { + /* Give each thread a disjoint, contiguous range of frequencies [klo,khi) (so the + * f_hat[k] writes are race-free) and run the phase recurrence within it, re-seeded + * every B steps. Cap the team at N/B threads so every thread owns at least one full + * block (>= B frequencies): with sub-block ranges the per-block seed (a COS+SIN pair) + * is un-amortised and the recurrence loses to the plain per-k loop. The cap depends + * only on N and B, so it removes any dependence on the machine's thread count. */ + const int nt_cap = (int)(ths->N_total / B); + #pragma omp parallel default(shared) num_threads(nt_cap) + { + const int nt = omp_get_num_threads(); + const int tid = omp_get_thread_num(); + const INT klo = (INT)(((long long)ths->N_total * tid) / nt); + const INT khi = (INT)(((long long)ths->N_total * (tid + 1)) / nt); + INT j; + for (j = 0; j < ths->M_total; j++) + { + const R x = ths->x[j]; + const R dphi = K2PI * x; + const C dw = COS(dphi) + II * SIN(dphi); + INT k_L = klo; + while (k_L < khi) + { + const R omega = X(reduced_omega)((R)(k_L + OFFSET), x); + C w = COS(omega) + II * SIN(omega); + INT kend = k_L + B; if (kend > khi) kend = khi; + for (; k_L < kend; k_L++) + { + f_hat[k_L] += f[j] * BASEPART(w); + w *= dw; + } + } + } + } + } + else + { + /* N <= B: the recurrence spans at most one block per thread-range, so its per-block + * seed/setup costs more than it saves once threaded. Use the plain per-k + * parallelisation. At these tiny N the per-entry phase error is in check. */ INT k_L; #pragma omp parallel for default(shared) private(k_L) for (k_L = 0; k_L < ths->N_total; k_L++) @@ -168,15 +239,43 @@ void X(adjoint_direct)(const X(plan) *ths) f_hat[k_L] += f[j] * BASE(omega); } } + } #else INT j; - for (j = 0; j < ths->M_total; j++) + if (ths->N_total > B) { - INT k_L; - for (k_L = 0; k_L < ths->N_total; k_L++) + for (j = 0; j < ths->M_total; j++) { - R omega = K2PI * ((R)(k_L + OFFSET)) * ths->x[j]; - f_hat[k_L] += f[j] * BASE(omega); + const R x = ths->x[j]; + const R dphi = K2PI * x; + const C dw = COS(dphi) + II * SIN(dphi); + INT k_L = 0; + while (k_L < ths->N_total) + { + const R omega = X(reduced_omega)((R)(k_L + OFFSET), x); + C w = COS(omega) + II * SIN(omega); + INT kend = k_L + B; if (kend > ths->N_total) kend = ths->N_total; + for (; k_L < kend; k_L++) + { + f_hat[k_L] += f[j] * BASEPART(w); + w *= dw; + } + } + } + } + else + { + /* N <= B: the recurrence spans at most one block, so its per-j seed/setup outweighs + * the saved transcendentals (which the plain loop vectorises); the un-reduced + * argument's phase error is still in check at these tiny N. */ + for (j = 0; j < ths->M_total; j++) + { + INT k_L; + for (k_L = 0; k_L < ths->N_total; k_L++) + { + R omega = K2PI * ((R)(k_L + OFFSET)) * ths->x[j]; + f_hat[k_L] += f[j] * BASE(omega); + } } } #endif diff --git a/tests/check.c b/tests/check.c index 0fc12551..38dd39d8 100644 --- a/tests/check.c +++ b/tests/check.c @@ -80,11 +80,14 @@ int main(void) #ifdef HAVE_NFCT #undef X #define X(name) NFCT(name) -#ifndef _OPENMP nfct = CU_add_suite("nfct", 0, 0); + /* The 1D direct/adjoint transforms are OpenMP-parallelised, so validate them in both the + serial and the threaded (checkall_threads) build. The remaining NFCT paths (fast, 2D/3D, + online) are exercised single-threaded only. */ CU_add_test(nfct, "nfct_1d_direct_file", X(check_1d_direct_file)); - CU_add_test(nfct, "nfct_1d_fast_file", X(check_1d_fast_file)); CU_add_test(nfct, "nfct_adjoint_1d_direct_file", X(check_adjoint_1d_direct_file)); +#ifndef _OPENMP + CU_add_test(nfct, "nfct_1d_fast_file", X(check_1d_fast_file)); CU_add_test(nfct, "nfct_adjoint_1d_fast_file", X(check_adjoint_1d_fast_file)); CU_add_test(nfct, "nfct_1d_online", X(check_1d_online)); CU_add_test(nfct, "nfct_adjoint_1d_online", X(check_adjoint_1d_online)); @@ -112,11 +115,13 @@ int main(void) #ifdef HAVE_NFST #undef X #define X(name) NFST(name) -#ifndef _OPENMP nfst = CU_add_suite("nfst", 0, 0); + /* As for NFCT: the 1D direct/adjoint transforms are OpenMP-parallelised, so validate them in + both serial and threaded builds; the fast, 2D/3D and online NFST paths stay serial-only. */ CU_add_test(nfst, "nfst_1d_direct_file", X(check_1d_direct_file)); - CU_add_test(nfst, "nfst_1d_fast_file", X(check_1d_fast_file)); CU_add_test(nfst, "nfst_adjoint_1d_direct_file", X(check_adjoint_1d_direct_file)); +#ifndef _OPENMP + CU_add_test(nfst, "nfst_1d_fast_file", X(check_1d_fast_file)); CU_add_test(nfst, "nfst_adjoint_1d_fast_file", X(check_adjoint_1d_fast_file)); CU_add_test(nfst, "nfst_1d_online", X(check_1d_online)); CU_add_test(nfst, "nfst_adjoint_1d_online", X(check_adjoint_1d_online));