From ddb90bf7abecb55f30ee60a5419a303c784e2585 Mon Sep 17 00:00:00 2001 From: Claude Date: Thu, 11 Jun 2026 00:10:40 +0000 Subject: [PATCH] perf(signal): opt-in FFT operator for the CIR ISTA solver (8-14x measured) MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Phi is a sub-DFT, so each ISTA mat-vec can run as one length-G FFT (O(G log G)) instead of a dense O(K*G) product — the dominant-latency-hazard finding from the beyond-SOTA optimization roadmap. New CirConfig::fft_operator, default FALSE: the dense path stays the bit-exact witness default. The FFT evaluates the same sums in a different order, so enabling it shifts float results in the last bits and requires regenerating any pinned witness — strictly opt-in per deployment. FftOperator (rustfft, planned once at CirEstimator::new, scratch buffers reused across the ISTA loop) dispatches inside ista_solve: Phi x = scale * forward-FFT(x) sampled at bins (k_idx mod G) Phi^H v = scale * unnormalised inverse-FFT of v scattered into those bins Warm-start and Lipschitz estimation stay dense at construction. Measured (criterion, same run, same machine): ht20: 2.22 ms -> 265 us (8.4x) ht40: 10.26 ms -> 717 us (14.3x) The real HE40 grid (K=484, G=1452) scales further per the O(K*G)/O(G log G) ratio. 3 new tests: FFT<->dense matvec equivalence to float tolerance on ht20 and he40 grids; end-to-end dominant-tap agreement on a single-path frame; all default configs keep FFT off. New cir_estimate_fft bench group. Workspace gate: 2,921 passed / 0 failed (default path bit-exact, witnesses unchanged). https://claude.ai/code/session_01MjBucx95K4BuUxZi8NWwRH --- CHANGELOG.md | 1 + .../benches/cir_bench.rs | 31 +++ .../wifi-densepose-signal/src/ruvsense/cir.rs | 210 +++++++++++++++++- 3 files changed, 239 insertions(+), 3 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 84bdeeb0..83a18207 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -8,6 +8,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ## [Unreleased] ### Added +- **Opt-in FFT operator for the CIR ISTA solver (8–14× measured).** Φ is a sub-DFT, so each ISTA mat-vec can run as one length-G FFT (O(G log G)) instead of a dense O(K·G) product. New `CirConfig::fft_operator` (default **false** — the dense path stays the bit-exact witness default; the FFT evaluates the same sums in a different order, so enabling it shifts float results and requires regenerating any pinned witness). `FftOperator` (rustfft, planned once at construction, scratch reused across the ISTA loop) dispatches inside `ista_solve`; warm-start/Lipschitz stay dense at construction. Measured (criterion, same run): ht20 2.22 ms → 265 µs (**8.4×**), ht40 10.26 ms → 717 µs (**14.3×**); the real HE40 grid (K=484, G=1452) scales further. 3 new tests: FFT↔dense matvec equivalence to float tolerance (ht20 + he40 grids), end-to-end dominant-tap agreement on a single-path frame, and all default configs keep FFT off. New `cir_estimate_fft` bench group. - **Per-room adapter provenance + drift→recalibration advisor in the streaming engine.** Closes the trust-chain gap where an ~11 KB per-room LoRA adapter (ADR-150 §3.4) could silently change inference without the witness noticing. `StreamingEngine::set_room_adapter(AdapterInfo)` pins the adapter's content-derived id into provenance `model_version` (`rfenc-v1+adapter:`) — and therefore into the BLAKE3 witness — so swapping or clearing adapter weights always shifts the witness (engine test proves base → adapter → other-adapter → cleared all witness differently, and cleared == base). New `RecalibrationAdvisor` recommends re-running the ADR-135 baseline / refitting the adapter on sustained low fusion coherence (streak threshold, default 60 cycles ≈ 3 s at 20 Hz) or an ADR-142 change-point; surfaced as `TrustedOutput::recalibration_recommended` and stored on the sensing-server `AppState` alongside the witness. Bridge plumbing: `EngineBridge::{set_room_adapter, clear_room_adapter}` + live-path test that the adapter id flows into the live witness. Engine 15 tests, bridge 7 tests. *Scope note: this is the deployable provenance/trigger half of the "retrained model" roadmap item — fitting the adapter itself runs in the existing external calibration service (`aether-arena/calibration/`), and a trained RF-encoder checkpoint still does not exist in-tree.* - **RuView beyond-SOTA research series** (`docs/research/ruview-beyond-sota/`, 6 docs) — research-swarm output defining the beyond-SOTA bar and the path to it: system capability audit (role→crate maturity matrix, gap analysis, risk register), web-verified 2026 SOTA landscape per capability axis (incl. ratified IEEE 802.11bf-2025), 8-pillar target architecture on the ADR-136 contract spine (no rewrite), 6-layer benchmark/validation methodology (all 15 criterion bench targets inventoried; ADR-149 statistical protocol), and a determinism-safe optimization roadmap. Includes session validation evidence: 2,797 workspace tests / 0 failed, Python proof PASS (bit-exact), paired pre/post criterion runs. diff --git a/v2/crates/wifi-densepose-signal/benches/cir_bench.rs b/v2/crates/wifi-densepose-signal/benches/cir_bench.rs index 1abf726e..30619759 100644 --- a/v2/crates/wifi-densepose-signal/benches/cir_bench.rs +++ b/v2/crates/wifi-densepose-signal/benches/cir_bench.rs @@ -156,6 +156,36 @@ fn bench_estimate(c: &mut Criterion) { group.finish(); } +// --------------------------------------------------------------------------- +// Benchmark 1b: opt-in FFT operator (CirConfig::fft_operator = true) +// --------------------------------------------------------------------------- + +/// Same workload as `cir_estimate`, with the O(G log G) FFT Φ/Φᴴ operator +/// enabled. Compare against `cir_estimate/` for the dense baseline. +fn bench_estimate_fft(c: &mut Criterion) { + let mut group = c.benchmark_group("cir_estimate_fft"); + + let tiers: &[(&str, u16)] = &[("ht20", 20), ("ht40", 40), ("he40", 40)]; + + for &(label, bw_mhz) in tiers { + let mut cfg = CirConfig::for_bandwidth_mhz(bw_mhz); + cfg.fft_operator = true; + let k_active = cfg.delay_bins / 3; + + group.throughput(Throughput::Elements(k_active as u64)); + + let est = CirEstimator::new(cfg.clone()); + let csi = synth_csi(&cfg); + let frame = make_frame(bw_mhz, csi); + + group.bench_with_input(BenchmarkId::from_parameter(label), &frame, |b, f| { + b.iter(|| black_box(est.estimate(black_box(f)).ok())); + }); + } + + group.finish(); +} + // --------------------------------------------------------------------------- // Benchmark 2: 12-link amortisation (shared estimator across links) // --------------------------------------------------------------------------- @@ -241,6 +271,7 @@ fn bench_estimator_construction(c: &mut Criterion) { criterion_group!( benches, bench_estimate, + bench_estimate_fft, bench_estimate_12link, bench_estimator_construction, ); diff --git a/v2/crates/wifi-densepose-signal/src/ruvsense/cir.rs b/v2/crates/wifi-densepose-signal/src/ruvsense/cir.rs index 0c692f50..d8669ff1 100644 --- a/v2/crates/wifi-densepose-signal/src/ruvsense/cir.rs +++ b/v2/crates/wifi-densepose-signal/src/ruvsense/cir.rs @@ -26,6 +26,8 @@ use num_complex::Complex32; use ruvector_solver::{neumann::NeumannSolver, types::CsrMatrix}; +use rustfft::{Fft, FftPlanner}; +use std::sync::Arc; use thiserror::Error; use wifi_densepose_core::types::CsiFrame; @@ -157,6 +159,16 @@ pub struct CirConfig { pub ranging_min_bw_hz: f64, /// Minimum dominant-tap ratio below which `ranging_valid` is false. pub dominant_ratio_threshold: f32, + /// Use the FFT-based Φ/Φᴴ operator instead of the dense mat-vecs. + /// + /// **Default `false` (dense, bit-exact witness path).** Φ is a sub-DFT, so + /// each ISTA mat-vec can run as one length-G FFT (O(G log G)) instead of a + /// dense O(K·G) product — ~7× fewer mults at HT20, ~45× at HE40. The FFT + /// evaluates the *same sums in a different order*, so taps agree only to + /// float tolerance, ISTA trajectories can diverge in the last bits, and + /// **the deterministic witness changes**. Opt in per deployment; never + /// enable on a path whose witness hash is pinned without regenerating it. + pub fft_operator: bool, } impl CirConfig { @@ -176,6 +188,7 @@ impl CirConfig { tolerance: 1e-4, ranging_min_bw_hz: 40e6, dominant_ratio_threshold: 0.3, + fft_operator: false, } } @@ -193,6 +206,7 @@ impl CirConfig { tolerance: 1e-4, ranging_min_bw_hz: 40e6, dominant_ratio_threshold: 0.3, + fft_operator: false, } } @@ -212,6 +226,7 @@ impl CirConfig { tolerance: 1e-4, ranging_min_bw_hz: 40e6, dominant_ratio_threshold: 0.3, + fft_operator: false, } } @@ -229,6 +244,7 @@ impl CirConfig { tolerance: 1e-4, ranging_min_bw_hz: 40e6, dominant_ratio_threshold: 0.3, + fft_operator: false, } } @@ -355,6 +371,87 @@ pub struct CirEstimator { warm_diag: Vec, /// Diagonal CSR matrix over `warm_diag` for the NeumannSolver warm-start. warm_csr: CsrMatrix, + /// FFT operator for Φ/Φᴴ, built only when `config.fft_operator` (opt-in). + fft: Option, +} + +/// FFT realisation of the sub-DFT sensing operator (opt-in, see +/// [`CirConfig::fft_operator`]). +/// +/// Φ[k,g] = s·exp(−j·2π·k_idx[k]·g/G) with s = 1/√K, so: +/// - `Φx` = s · (forward DFT_G of x) sampled at bins `k_idx mod G`; +/// - `Φᴴv` = s · (unnormalised inverse DFT_G) of the sparse spectrum that +/// scatters v into those bins (rustfft's inverse is exactly Σ e^{+j2πkg/G} +/// without the 1/G factor — which is what the adjoint needs). +/// +/// Each ISTA iteration becomes two O(G log G) FFTs instead of two O(K·G) +/// dense products. +struct FftOperator { + forward: Arc>, + inverse: Arc>, + /// Active-subcarrier DFT bins: `k_idx mod G`, one per active subcarrier. + bins: Vec, + /// 1/√K column normalisation of Φ. + scale: f32, + g: usize, +} + +impl FftOperator { + fn new(active_indices: &[i32], g: usize, k: usize) -> Self { + let mut planner = FftPlanner::::new(); + let bins = active_indices + .iter() + .map(|&idx| (idx.rem_euclid(g as i32)) as usize) + .collect(); + Self { + forward: planner.plan_fft_forward(g), + inverse: planner.plan_fft_inverse(g), + bins, + scale: 1.0 / (k as f32).sqrt(), + g, + } + } + + /// Φ v → out (out length K). `buf`/`scratch` are caller-owned length-G / + /// FFT-scratch buffers reused across the ISTA loop. + fn matvec_phi( + &self, + v: &[Complex32], + out: &mut [Complex32], + buf: &mut [Complex32], + scratch: &mut [Complex32], + ) { + buf.copy_from_slice(v); + self.forward.process_with_scratch(buf, scratch); + for (o, &bin) in out.iter_mut().zip(&self.bins) { + *o = buf[bin] * self.scale; + } + } + + /// Φᴴ v → out (out length G). + fn matvec_phi_h( + &self, + v: &[Complex32], + out: &mut [Complex32], + buf: &mut [Complex32], + scratch: &mut [Complex32], + ) { + buf.fill(Complex32::new(0.0, 0.0)); + for (&vi, &bin) in v.iter().zip(&self.bins) { + buf[bin] += vi; + } + self.inverse.process_with_scratch(buf, scratch); + for (o, &b) in out.iter_mut().zip(buf.iter()) { + *o = b * self.scale; + } + } + + /// Length of the FFT scratch buffer required by both plans. + fn scratch_len(&self) -> usize { + self.forward + .get_inplace_scratch_len() + .max(self.inverse.get_inplace_scratch_len()) + } } // Φ and Φ^H are immutable after construction; all `estimate()` locals are @@ -371,6 +468,9 @@ impl CirEstimator { let (phi, phi_h) = build_sensing_matrix(&active_indices, g, k); let lipschitz = estimate_lipschitz(&phi, &phi_h, k, g, 30); let (warm_diag, warm_csr) = build_warm_start_system(&phi, k, g, config.lambda); + let fft = config + .fft_operator + .then(|| FftOperator::new(&active_indices, g, k)); Self { config, sensing_matrix: phi, @@ -379,6 +479,7 @@ impl CirEstimator { lipschitz, warm_diag, warm_csr, + fft, } } @@ -420,6 +521,7 @@ impl CirEstimator { self.lipschitz, &self.warm_diag, &self.warm_csr, + self.fft.as_ref(), )?; let tap_sum: f32 = x.iter().map(|c| c.norm()).sum(); @@ -617,6 +719,7 @@ fn ista_solve( lipschitz: f32, warm_diag: &[f32], warm_csr: &CsrMatrix, + fft: Option<&FftOperator>, ) -> Result<(Vec, u32, f32), CirError> { let k = config.num_active; let g = config.num_taps; @@ -627,16 +730,31 @@ fn ista_solve( let mut x_prev = x.clone(); let mut phi_x = vec![Complex32::new(0.0, 0.0); k]; let mut grad = vec![Complex32::new(0.0, 0.0); g]; + // FFT-path work buffers, allocated once per solve (not per iteration). + let (mut fft_buf, mut fft_scratch) = match fft { + Some(op) => ( + vec![Complex32::new(0.0, 0.0); op.g], + vec![Complex32::new(0.0, 0.0); op.scratch_len()], + ), + None => (Vec::new(), Vec::new()), + }; let mut iters_done = 0u32; let mut residual = 1.0_f32; for iter in 0..config.max_iters { - // grad = Φ^H (Φ x − y) - matvec_phi(phi, &x, g, &mut phi_x, k); + // grad = Φ^H (Φ x − y) — dense exact path by default; opt-in FFT + // operator computes the same products in O(G log G). + match fft { + Some(op) => op.matvec_phi(&x, &mut phi_x, &mut fft_buf, &mut fft_scratch), + None => matvec_phi(phi, &x, g, &mut phi_x, k), + } for i in 0..k { phi_x[i] -= y[i]; } - matvec_phi_h(phi_h, &phi_x, k, &mut grad, g); + match fft { + Some(op) => op.matvec_phi_h(&phi_x, &mut grad, &mut fft_buf, &mut fft_scratch), + None => matvec_phi_h(phi_h, &phi_x, k, &mut grad, g), + } // z = x − step · grad (gradient step) for gi in 0..g { @@ -1049,4 +1167,90 @@ mod tests { let meta = CsiMetadata::new(DeviceId::new("test"), FrequencyBand::Band2_4GHz, 6); CsiFrame::new(meta, data) } + + // ---- Opt-in FFT operator (CirConfig::fft_operator) ---- + + /// The FFT operator computes the same Φ/Φᴴ products as the dense path to + /// float tolerance, for both a small (HT20) and the largest (HE40) config. + #[test] + fn fft_matvecs_match_dense() { + for config in [CirConfig::ht20(), CirConfig::he40()] { + let k = config.num_active; + let g = config.num_taps; + let active: Vec = config.active_indices().to_vec(); + let (phi, phi_h) = build_sensing_matrix(&active, g, k); + let op = FftOperator::new(&active, g, k); + let mut buf = vec![Complex32::new(0.0, 0.0); g]; + let mut scratch = vec![Complex32::new(0.0, 0.0); op.scratch_len()]; + + // Deterministic non-trivial input vectors. + let x: Vec = (0..g) + .map(|i| Complex32::new((i as f32 * 0.37).sin(), (i as f32 * 0.71).cos())) + .collect(); + let v: Vec = (0..k) + .map(|i| Complex32::new((i as f32 * 0.13).cos(), (i as f32 * 0.29).sin())) + .collect(); + + // Φx: dense vs FFT. + let mut dense_kx = vec![Complex32::new(0.0, 0.0); k]; + matvec_phi(&phi, &x, g, &mut dense_kx, k); + let mut fft_kx = vec![Complex32::new(0.0, 0.0); k]; + op.matvec_phi(&x, &mut fft_kx, &mut buf, &mut scratch); + let scale_ref: f32 = dense_kx.iter().map(|c| c.norm()).sum::() / k as f32; + for (d, f) in dense_kx.iter().zip(&fft_kx) { + assert!( + (d - f).norm() <= 1e-3 * scale_ref.max(1.0), + "phi matvec mismatch (G={g}): {d} vs {f}" + ); + } + + // Φᴴv: dense vs FFT. + let mut dense_gv = vec![Complex32::new(0.0, 0.0); g]; + matvec_phi_h(&phi_h, &v, k, &mut dense_gv, g); + let mut fft_gv = vec![Complex32::new(0.0, 0.0); g]; + op.matvec_phi_h(&v, &mut fft_gv, &mut buf, &mut scratch); + let scale_ref_g: f32 = dense_gv.iter().map(|c| c.norm()).sum::() / g as f32; + for (d, f) in dense_gv.iter().zip(&fft_gv) { + assert!( + (d - f).norm() <= 1e-3 * scale_ref_g.max(1.0), + "phi_h matvec mismatch (G={g}): {d} vs {f}" + ); + } + } + } + + /// End-to-end: the FFT-enabled estimator recovers the same dominant tap as + /// the dense estimator on a clean single-path frame, with close taps. + #[test] + fn fft_estimate_matches_dense_dominant_tap() { + let dense_cfg = CirConfig::ht20(); + let mut fft_cfg = CirConfig::ht20(); + fft_cfg.fft_operator = true; + + let frame = make_single_tap_frame(dense_cfg.num_subcarriers, 50e-9); + let dense = CirEstimator::new(dense_cfg).estimate(&frame).unwrap(); + let fast = CirEstimator::new(fft_cfg).estimate(&frame).unwrap(); + + assert_eq!(dense.dominant_tap_idx, fast.dominant_tap_idx); + assert!((dense.dominant_tap_ratio - fast.dominant_tap_ratio).abs() < 1e-2); + // Tap vectors agree to float tolerance relative to the dominant tap. + let dom = dense.taps[dense.dominant_tap_idx].norm().max(1e-6); + for (a, b) in dense.taps.iter().zip(&fast.taps) { + assert!((a - b).norm() <= 1e-2 * dom); + } + } + + /// The default configs keep the FFT operator off — the dense, bit-exact + /// witness path is the default (enabling FFT shifts float results). + #[test] + fn fft_operator_is_off_by_default() { + for c in [ + CirConfig::ht20(), + CirConfig::ht40(), + CirConfig::he20(), + CirConfig::he40(), + ] { + assert!(!c.fft_operator); + } + } }