|
| 1 | +#!/usr/bin/env python3 |
| 2 | +""" |
| 3 | +LowPassFilter2p correspondence tests — Route B. |
| 4 | +
|
| 5 | +🔬 Lean Squad automated correspondence validation. |
| 6 | +
|
| 7 | +Validates that the Lean 4 model in FVSquad/LowPassFilter2p.lean |
| 8 | +matches the C++ implementation in src/lib/mathlib/math/filter/LowPassFilter2p.hpp |
| 9 | +on a shared set of test cases. |
| 10 | +
|
| 11 | +C++ logic (Direct Form II): |
| 12 | + T apply(const T &sample) { |
| 13 | + T delay_element_0 = sample - _delay_element_1 * _a1 - _delay_element_2 * _a2; |
| 14 | + T output = delay_element_0 * _b0 + _delay_element_1 * _b1 + _delay_element_2 * _b2; |
| 15 | + _delay_element_2 = _delay_element_1; |
| 16 | + _delay_element_1 = delay_element_0; |
| 17 | + return output; |
| 18 | + } |
| 19 | +
|
| 20 | +Lean model (over Rat/fractions): |
| 21 | + def lpf2pApply (c : LPF2pCoeffs) (s : LPF2pState) (sample : Rat) : LPF2pState × Rat := |
| 22 | + let w0 := sample - s.d1 * c.a1 - s.d2 * c.a2 |
| 23 | + let out := w0 * c.b0 + s.d1 * c.b1 + s.d2 * c.b2 |
| 24 | + ({ d1 := w0, d2 := s.d1 }, out) |
| 25 | +
|
| 26 | +Since the Lean model uses exact rational arithmetic and the C++ uses float, |
| 27 | +we compare using Python fractions (exact rational arithmetic) to verify algorithmic |
| 28 | +correspondence, isolating model structure from floating-point artefacts. |
| 29 | +
|
| 30 | +Run with: python3 test_correspondence.py |
| 31 | +""" |
| 32 | + |
| 33 | +from fractions import Fraction |
| 34 | +import sys |
| 35 | +from dataclasses import dataclass |
| 36 | +from typing import List, Tuple |
| 37 | + |
| 38 | + |
| 39 | +# ── Lean model implementation (Python/Fraction) ────────────────────────────── |
| 40 | + |
| 41 | +@dataclass |
| 42 | +class LPF2pState: |
| 43 | + d1: Fraction |
| 44 | + d2: Fraction |
| 45 | + |
| 46 | +@dataclass |
| 47 | +class LPF2pCoeffs: |
| 48 | + b0: Fraction |
| 49 | + b1: Fraction |
| 50 | + b2: Fraction |
| 51 | + a1: Fraction |
| 52 | + a2: Fraction |
| 53 | + |
| 54 | +def lean_lpf2p_apply(c: LPF2pCoeffs, s: LPF2pState, sample: Fraction) -> Tuple[LPF2pState, Fraction]: |
| 55 | + """Direct translation of lpf2pApply from LowPassFilter2p.lean.""" |
| 56 | + w0 = sample - s.d1 * c.a1 - s.d2 * c.a2 |
| 57 | + out = w0 * c.b0 + s.d1 * c.b1 + s.d2 * c.b2 |
| 58 | + return LPF2pState(d1=w0, d2=s.d1), out |
| 59 | + |
| 60 | +def lean_apply_n(c: LPF2pCoeffs, s0: LPF2pState, samples: List[Fraction]) -> Tuple[LPF2pState, List[Fraction]]: |
| 61 | + """Apply lean model for a sequence of samples, returning final state and all outputs.""" |
| 62 | + s = s0 |
| 63 | + outputs = [] |
| 64 | + for sample in samples: |
| 65 | + s, out = lean_lpf2p_apply(c, s, sample) |
| 66 | + outputs.append(out) |
| 67 | + return s, outputs |
| 68 | + |
| 69 | + |
| 70 | +# ── C++ model implementation (Python/Fraction, same algorithm) ─────────────── |
| 71 | + |
| 72 | +def cpp_lpf2p_apply(c: LPF2pCoeffs, s: LPF2pState, sample: Fraction) -> Tuple[LPF2pState, Fraction]: |
| 73 | + """Direct translation of C++ apply() — using Fraction for exact comparison.""" |
| 74 | + delay_element_0 = sample - s.d1 * c.a1 - s.d2 * c.a2 |
| 75 | + output = delay_element_0 * c.b0 + s.d1 * c.b1 + s.d2 * c.b2 |
| 76 | + new_d2 = s.d1 |
| 77 | + new_d1 = delay_element_0 |
| 78 | + return LPF2pState(d1=new_d1, d2=new_d2), output |
| 79 | + |
| 80 | +def cpp_apply_n(c: LPF2pCoeffs, s0: LPF2pState, samples: List[Fraction]) -> Tuple[LPF2pState, List[Fraction]]: |
| 81 | + s = s0 |
| 82 | + outputs = [] |
| 83 | + for sample in samples: |
| 84 | + s, out = cpp_lpf2p_apply(c, s, sample) |
| 85 | + outputs.append(out) |
| 86 | + return s, outputs |
| 87 | + |
| 88 | + |
| 89 | +# ── Test cases ──────────────────────────────────────────────────────────────── |
| 90 | + |
| 91 | +def F(x): |
| 92 | + """Convert string or number to Fraction.""" |
| 93 | + return Fraction(x) |
| 94 | + |
| 95 | +def make_coeffs(b0, b1, b2, a1, a2): |
| 96 | + return LPF2pCoeffs(F(b0), F(b1), F(b2), F(a1), F(a2)) |
| 97 | + |
| 98 | +def zero_state(): |
| 99 | + return LPF2pState(F(0), F(0)) |
| 100 | + |
| 101 | +# Disabled filter: b0=1, b1=b2=a1=a2=0 |
| 102 | +DISABLED = make_coeffs(1, 0, 0, 0, 0) |
| 103 | + |
| 104 | +# Butterworth-like rational approximation (not a real Butterworth, but exercises b1=2*b0, b2=b0) |
| 105 | +# Using b0=1/4, b1=1/2, b2=1/4, a1=-1/2, a2=1/4 (arbitrary rational coefficients) |
| 106 | +BUTTERWORTH_RAT = make_coeffs("1/4", "1/2", "1/4", "-1/2", "1/4") |
| 107 | + |
| 108 | +# Another set: b0=1/10, b1=2/10, b2=1/10, a1=-3/5, a2=1/5 |
| 109 | +COEFF_B = make_coeffs("1/10", "2/10", "1/10", "-3/5", "1/5") |
| 110 | + |
| 111 | +# DC-unity coefficients: b0+b1+b2 = 1+a1+a2 |
| 112 | +# b0=1/6, b1=2/6, b2=1/6, a1=-1/6, a2=-1/6: sum(b)=4/6, sum(a_den)=1+(-1/6)+(-1/6)=4/6 ✓ |
| 113 | +DC_UNITY = make_coeffs("1/6", "2/6", "1/6", "-1/6", "-1/6") |
| 114 | + |
| 115 | +CASES = [ |
| 116 | + # (name, coeffs, d1_init, d2_init, samples) |
| 117 | + # Basic disabled pass-through |
| 118 | + ("disabled_zero_state_single", DISABLED, "0", "0", ["1"]), |
| 119 | + ("disabled_zero_state_multi", DISABLED, "0", "0", ["1", "2", "3", "-1", "0"]), |
| 120 | + ("disabled_nonzero_state", DISABLED, "5", "-3", ["7", "-2", "0"]), |
| 121 | + ("disabled_neg_input", DISABLED, "0", "0", ["-5", "-10", "-15"]), |
| 122 | + ("disabled_zero_input", DISABLED, "3", "7", ["0", "0", "0"]), |
| 123 | + |
| 124 | + # Butterworth-like rational coefficients |
| 125 | + ("bw_rat_zero_in", BUTTERWORTH_RAT, "0", "0", ["0"]*5), |
| 126 | + ("bw_rat_unit_step", BUTTERWORTH_RAT, "0", "0", ["1"]*8), |
| 127 | + ("bw_rat_neg_step", BUTTERWORTH_RAT, "0", "0", ["-1"]*8), |
| 128 | + ("bw_rat_alternating", BUTTERWORTH_RAT, "0", "0", ["1", "-1"]*5), |
| 129 | + ("bw_rat_nonzero_init", BUTTERWORTH_RAT, "1/2", "-1/3", ["1", "2", "3"]), |
| 130 | + ("bw_rat_large_input", BUTTERWORTH_RAT, "0", "0", ["100", "100", "100"]), |
| 131 | + ("bw_rat_fractional_input", BUTTERWORTH_RAT, "0", "0", ["1/3", "2/3", "1"]), |
| 132 | + |
| 133 | + # Coefficient set B |
| 134 | + ("coeff_b_zero_in", COEFF_B, "0", "0", ["0"]*5), |
| 135 | + ("coeff_b_unit_step", COEFF_B, "0", "0", ["1"]*10), |
| 136 | + ("coeff_b_impulse", COEFF_B, "0", "0", ["1"] + ["0"]*9), |
| 137 | + ("coeff_b_ramp", COEFF_B, "0", "0", [str(k) for k in range(10)]), |
| 138 | + ("coeff_b_nonzero_init", COEFF_B, "3/4", "1/4", ["0"]*5), |
| 139 | + |
| 140 | + # DC unity: DC gain = 1 → steady-state output should equal input |
| 141 | + ("dc_unity_zero_state_const", DC_UNITY, "0", "0", ["1"]*20), |
| 142 | + ("dc_unity_zero_state_neg", DC_UNITY, "0", "0", ["-1"]*20), |
| 143 | + ("dc_unity_impulse", DC_UNITY, "0", "0", ["1"] + ["0"]*15), |
| 144 | + ("dc_unity_step_various", DC_UNITY, "0", "0", ["1/2"]*15), |
| 145 | + |
| 146 | + # Zero input, nonzero state (exponential decay) |
| 147 | + ("bw_rat_zero_in_nonzero_state", BUTTERWORTH_RAT, "1", "1", ["0"]*10), |
| 148 | + ("coeff_b_zero_in_nonzero_state", COEFF_B, "2", "-1", ["0"]*10), |
| 149 | + |
| 150 | + # Single-sample checks (boundary) |
| 151 | + ("disabled_single_pos", DISABLED, "0", "0", ["42"]), |
| 152 | + ("disabled_single_neg", DISABLED, "0", "0", ["-7"]), |
| 153 | + ("disabled_single_frac", DISABLED, "0", "0", ["1/7"]), |
| 154 | + ("bw_rat_single", BUTTERWORTH_RAT, "0", "0", ["5"]), |
| 155 | + ("coeff_b_single", COEFF_B, "1/2", "1/3", ["2/5"]), |
| 156 | +] |
| 157 | + |
| 158 | + |
| 159 | +# ── Test runner ─────────────────────────────────────────────────────────────── |
| 160 | + |
| 161 | +def run_tests(): |
| 162 | + passed = 0 |
| 163 | + failed = 0 |
| 164 | + errors = [] |
| 165 | + |
| 166 | + for case in CASES: |
| 167 | + name, c, d1_str, d2_str, sample_strs = case |
| 168 | + d1 = F(d1_str) |
| 169 | + d2 = F(d2_str) |
| 170 | + samples = [F(s) for s in sample_strs] |
| 171 | + |
| 172 | + s0 = LPF2pState(d1=d1, d2=d2) |
| 173 | + |
| 174 | + lean_s, lean_outs = lean_apply_n(c, s0, samples) |
| 175 | + cpp_s, cpp_outs = cpp_apply_n(c, s0, samples) |
| 176 | + |
| 177 | + ok = (lean_outs == cpp_outs and |
| 178 | + lean_s.d1 == cpp_s.d1 and |
| 179 | + lean_s.d2 == cpp_s.d2) |
| 180 | + |
| 181 | + if ok: |
| 182 | + passed += 1 |
| 183 | + else: |
| 184 | + failed += 1 |
| 185 | + diff_outs = [(i, l, cr) for i, (l, cr) in enumerate(zip(lean_outs, cpp_outs)) if l != cr] |
| 186 | + errors.append((name, diff_outs, lean_s, cpp_s)) |
| 187 | + |
| 188 | + if errors: |
| 189 | + print(f"FAILED: {failed}/{passed + failed} test cases failed\n") |
| 190 | + for name, diff_outs, lean_s, cpp_s in errors: |
| 191 | + print(f" CASE: {name}") |
| 192 | + for i, l, cr in diff_outs: |
| 193 | + print(f" output[{i}]: lean={l}, cpp={cr}") |
| 194 | + if lean_s.d1 != cpp_s.d1: |
| 195 | + print(f" final d1: lean={lean_s.d1}, cpp={cpp_s.d1}") |
| 196 | + if lean_s.d2 != cpp_s.d2: |
| 197 | + print(f" final d2: lean={lean_s.d2}, cpp={cpp_s.d2}") |
| 198 | + sys.exit(1) |
| 199 | + else: |
| 200 | + total = passed + failed |
| 201 | + print(f"OK: {total}/{total} test cases passed") |
| 202 | + print(f"Tested: disabled pass-through, Butterworth-like rational, DC-unity, impulse response, ramp, zero-input decay") |
| 203 | + sys.exit(0) |
| 204 | + |
| 205 | + |
| 206 | +if __name__ == "__main__": |
| 207 | + run_tests() |
0 commit comments