Skip to content

Add OnePoleBank: phasor-free resonator core (~1.4–1.7× faster)#6

Open
enomado wants to merge 3 commits into
jhartquist:mainfrom
enomado:onepole-bank
Open

Add OnePoleBank: phasor-free resonator core (~1.4–1.7× faster)#6
enomado wants to merge 3 commits into
jhartquist:mainfrom
enomado:onepole-bank

Conversation

@enomado

@enomado enomado commented Jun 10, 2026

Copy link
Copy Markdown

What

Adds OnePoleBank, a reformulation of the resonator core that drops the per-sample running phasor and the periodic stabilization it requires, while producing bit-identical power/magnitude. It's ~1.35–1.70× faster than the current ResonatorBank on x86-64-v3 (AVX2), and is numerically self-stabilizing.

It's added alongside ResonatorBank — purely additive, no behavior change to existing code. Happy to fold it into ResonatorBank or make it the default core if you like the approach.

Why it works

The current resonator keeps, per bin, a unit phasor z rotated by e^{-iωΔt} every sample, accumulates α·x·z into the EWMA, and must periodically renormalize z because |w| != 1 exactly in f32 (so z drifts) — hence STABILIZE_EVERY and the sqrt renorm.

A resonator is a heterodyne followed by a one-pole lowpass. Move into the de-rotated frame u_n = r_n · e^{+iωnΔt} and the "rotate phasor + EWMA" pair collapses into a single complex one-pole with a constant complex coefficient:

u_n = c·u_{n-1} + α·x_n      c = (1-α)·e^{+iωΔt}
v_n = d·v_{n-1} + β·u_n      d = (1-β)·e^{+iωΔt}   (the output-smoothing EWMA)

Two facts make this a strict win:

  1. Power/magnitude are unchanged. |e^{-iωnΔt}| = 1, so |v_n| ≡ |rr_n| exactly. For spectrogram/magnitude use the phasor was never needed.
  2. It's self-stabilizing. |c| = (1-α) < 1 strictly, so the recurrence is contractive — round-off decays instead of accumulating. There is nothing to renormalize. STABILIZE_EVERY, the sqrt, and the chunk-around-stabilize loop in process_samples all go away.

The absolute-frame complex value (and phase) is still available — recovered by de-rotating v by e^{-iωnΔt} at readout only (hop boundaries), instead of maintaining a phasor every sample. The de-rotation argument grows without bound, so it's reduced in f64 (see Correctness).

Correctness (all in this PR)

Power/magnitude is algebraically equal to ResonatorBank for any (freq, α, β) — only f32 round-off separates them — so the bar is parity, not "close enough":

  • Fixture paritytests/reference_onepole.rs runs against the same chirp_88.npz fixture as reference.rs:
    • power matches to 9.9e-5 relative (pure f32 round-off — the two are mathematically identical),
    • complex matches to 6.1e-4 via lazy de-rotation.
  • Property-based parity — a proptest fuzzes the consumer's operating point: independent alpha_scale/beta_scale (so α != β), the full MIDI 12..=84 range including the high-Q low end where heuristic_alpha is tiny, and an arbitrary multi-partial stimulus (96 cases), asserting OnePoleBank powers track ResonatorBank within 0.5 % rel + 0.1 % peak. Worst case measured on a real 5-bins/semitone grid: ~6e-5 relative (~1e-4 with extreme α != β) — far below any display threshold, so a ResonatorBank → OnePoleBank swap cannot change a magnitude/power readout.
  • Long-run phase stabilityphase()/complex() de-rotate v by e^{-iωnΔt}, whose argument grows without bound; it's formed and reduced in f64 so a long-running readout stays accurate (at f32 the trig argument degrades within minutes). phase_readout_stable_at_large_sample_count recovers a +1 Hz detuning from two phase samples 128 apart after ~4 M samples — which fails outright if the de-rotation is done in f32.
  • Inline test asserts power equality vs ResonatorBank on a multi-tone signal.

Full crate suite green: 27 unit (incl. the property test) + 2 + 2 fixture/integration + 2 doc.

SIMD preserved

Same SoA layout and the same mul_add (with the existing wasm32+simd128 unfuse). Confirmed in the release asm that the hot loop autovectorizes to AVX2 vfmadd231ps over ymm (8 resonators/iteration) reading the strided SoA arrays — i.e. it vectorizes the same way ResonatorBank does.

Benchmark

benches/compare.rs benches both banks side by side (1 s signal, per-sample). On x86-64-v3 / AVX2:

bins phasor (current) onepole speedup
88 29.0 M/s 49.2 M/s 1.70×
264 12.0 M/s 17.5 M/s 1.46×
440 7.4 M/s 10.6 M/s 1.44×
880 3.6 M/s 4.8 M/s 1.35×

The speedup is largest when the working set is cache-resident (compute-bound, small bin counts) and shrinks as it becomes memory-bandwidth-bound — consistent with the change being a per-sample arithmetic reduction (the phasor rotate is gone).

Files

  • src/bank_onepole.rsOnePoleBank (full read API + inline + property-based parity tests)
  • src/lib.rs — export
  • tests/reference_onepole.rs — fixture parity (power + complex)
  • benches/compare.rs + Cargo.toml — side-by-side bench
  • Cargo.tomlproptest dev-dependency (property-based parity)

🤖 Generated with Claude Code

enomado and others added 2 commits June 10, 2026 03:02
Reformulate the resonator update into two cascaded complex one-poles with
constant complex coefficients, eliminating the per-sample running phasor and
the periodic stabilization it requires.

The reference ResonatorBank maintains, per resonator, a unit phasor z that it
rotates by e^{-iωΔt} every sample, and must periodically renormalize because
|w| != 1 exactly in f32 (z drifts). Working in the de-rotated frame
u_n = r_n·e^{+iωnΔt} collapses "rotate phasor + EWMA" into a single complex
one-pole:

    u_n = c·u_{n-1} + α·x_n      c = (1-α)·e^{+iωΔt}
    v_n = d·v_{n-1} + β·u_n      d = (1-β)·e^{+iωΔt}   (output smoothing)

Because |e^{-iωnΔt}| = 1, power and magnitude are identical to the reference
(|v_n| ≡ |rr_n|). The absolute-frame complex value / phase is recovered by
de-rotating at readout only (hop boundaries), never per sample.

Numerically this is strictly contractive: |c| = (1-α) < 1, so round-off decays
instead of accumulating. There is nothing to stabilize — STABILIZE_EVERY, the
sqrt renormalization, and the chunk-around-stabilize logic all disappear.

Verification (in this commit):
- tests/reference_onepole.rs checks against the same chirp_88 fixture: power
  matches to 9.9e-5 relative (f32 round-off), complex to 6.1e-4 via de-rotation.
- Inline test asserts power equality vs ResonatorBank on a multi-tone signal.
- Hot loop still autovectorizes: AVX2 vfmadd231ps over ymm confirmed in asm.

benches/compare.rs (M2 numbers will vary; x86-64-v3 / AVX2):
    bins  phasor   onepole   speedup
     88   29.0 M/s 49.2 M/s   1.70×
    264   12.0 M/s 17.5 M/s   1.46×
    440    7.4 M/s 10.6 M/s   1.44×
    880    3.6 M/s  4.8 M/s   1.35×

OnePoleBank is added alongside ResonatorBank (additive, no behavior change). It
could become the default core, or fold into ResonatorBank, at your discretion.

Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
complex()/phase() recover the absolute-frame value by de-rotating the stored
state by e^{-iωnΔt}. The angle 2π·f·n/sr grows without bound as n accumulates;
at f32 the sin/cos argument degrades within minutes (ulp ≈ 0.2 rad after ~10 min
at 440 Hz) and is meaningless after ~1 hour, so long-running phase/complex
readouts drift into noise. Form and reduce the angle in f64 instead — accurate
for many hours. Power is frame-invariant and unaffected.

Adds a regression test: detuning recovered from two phase samples 128 apart
after ~90 s of audio (n ≈ 4M, trig argument ~2.5e6 rad) still matches the true
+1 Hz offset. This fails outright with f32 de-rotation.

Motivated by driving a live pitch-spiral tuner off phase(), where the f32
formulation visibly drifts over a session.

Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
The reformulation's core claim is |v_n| ≡ |rr_n|: power and magnitude are
algebraically identical to the phasor bank for any (freq, alpha, beta), and
the f64 de-rotation only affects phase()/complex(). The existing
`power_matches_reference_bank` pins that at a single operating point
(alpha == beta, 5 fixed freqs, one signal).

Promote it to a property: fuzz independent alpha_scale/beta_scale (so
alpha != beta), the full MIDI 12..=84 range including the high-Q low end
where heuristic_alpha is tiny, and an arbitrary multi-partial stimulus
(96 cases). Asserts powers agree within 0.5% rel + 0.1% peak.

Measured worst case on the consumer's real 5-bins/semitone grid is ~6e-5
relative (~1e-4 with extreme alpha!=beta) — far below any display
threshold, so a ResonatorBank -> OnePoleBank swap cannot change a
magnitude/power readout. Adds proptest as a dev-dependency.

Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant