From 02421ae379244134baaf0850130e54e297b3ad51 Mon Sep 17 00:00:00 2001 From: Claude Date: Fri, 30 Oct 2026 11:54:52 +0000 Subject: [PATCH] =?UTF-8?q?research:=20D=C3=BCtsch=20analyzer=20=E2=80=94?= =?UTF-8?q?=20explicit=20opt-in=20horizon,=20complete-cycle=20check,=20HF?= =?UTF-8?q?=20component=20report?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Default (no --horizon) keeps the original cap-11 scoring. --horizon N requires complete cycles 1..N in every case, passes the horizon to common_window and labels the JSON. The unqualified phase-locked HF force component is reported per case. Protocol section 6b declares the D40/60D 31-cycle duration-only diagnostic. Co-Authored-By: Claude Opus 5.5 --- docs/validation/duetsch_re100_kc5_protocol.md | 19 ++++++++ research/duetsch_re100_kc5/analyze.py | 43 ++++++++++++++++--- research/duetsch_re100_kc5/protocol.py | 4 +- tests/unit/test_duetsch_protocol.py | 4 ++ 4 files changed, 62 insertions(+), 8 deletions(-) diff --git a/docs/validation/duetsch_re100_kc5_protocol.md b/docs/validation/duetsch_re100_kc5_protocol.md index 8421d15..d4ff678 100644 --- a/docs/validation/duetsch_re100_kc5_protocol.md +++ b/docs/validation/duetsch_re100_kc5_protocol.md @@ -204,6 +204,25 @@ Superseded (earlier drafts): the velocity-jump start, guessed stations and phases, "C_M" for c_i, and the uncertainty wording; cycles 4–6, the 30-cycle cap and the earlier tolerance numbers. +## 6b. Amendment: duration-only diagnostic (declared before the run) + +- Background: the D40 domain stage (20 D / 40 D / 60 D, 11.25 T each) failed + settling in every box (no window). That result stands. +- One declared diagnostic: D = 40, 60 D box, from rest, through cycle 31 + (157,500 steps = 31.25 T). Physics, settling tolerances, bands and + production unchanged; three proof pairs plus a continuously settled + K..K+2 with K ≤ 29. +- Invocation: `run.py --D 40 --box 60 --cycles 31` then + `analyze.py --horizon 31`; the JSON records `horizon` and `horizon_mode`. + Without `--horizon` the analyzer is the original cap-11 scoring. +- Failure at 31.25 T: report and stop; no automatic extension, matrix, D80 + or VIV. A pass shows only that this box settles with time; it is not a + domain or literature pass. +- The phase-locked high-frequency force component (residual of the mean + waveform after its Morison reconstruction, above 20/T) is reported per case + as unqualified; it needs separate qualification before any Morison or + literature statement. + ## 6a. Driver unit tests (CPU, before any run) - Phase/cycle arithmetic: s ↔ (k, φ), cycle boundaries, bin membership, diff --git a/research/duetsch_re100_kc5/analyze.py b/research/duetsch_re100_kc5/analyze.py index 10d8ec8..0061fbc 100644 --- a/research/duetsch_re100_kc5/analyze.py +++ b/research/duetsch_re100_kc5/analyze.py @@ -92,8 +92,28 @@ def pair_metrics(case, k): return m -def ledger(case): - ks = sorted(cycle_slices(case)) +def hf_component(case, ks): + """Unqualified phase-locked high-frequency force component (reported, not checked): + mean CME waveform over cycles ks minus its two-basis Morison reconstruction, + split at 20/T; cycle differencing does not remove it.""" + T = case["T"] + mean = np.mean([coarse_wave(case, k) for k in ks], axis=0) + c_d, c_i = cycle_morison(case, ks) + q = T // 4 + b = np.array([p.interval_bases(q + m, T) for m in range(T)]) + resid = mean - (c_d * b[:, 0] - (np.pi ** 2 * c_i / p.KC) * b[:, 1]) + spec = np.fft.rfft(resid) + high = spec.copy() + high[:21] = 0.0 + energy = np.abs(spec) ** 2 + return dict(label="unqualified phase-locked HF force component, not checked", + cycles=list(ks), hf_rms_over_F_ref=float(np.sqrt(np.mean(np.fft.irfft(high, T) ** 2))), + hf_energy_fraction=float(energy[21:].sum() / max(energy.sum(), 1e-300)), + top_harmonics=[int(j) for j in np.argsort(energy[21:])[::-1][:5] + 21]) + + +def ledger(case, horizon=p.LAST_SCORED): + ks = [k for k in sorted(cycle_slices(case)) if k <= horizon] pairs = {k: pair_metrics(case, k) for k in ks if k - 1 > p.DISCARD} return pairs, p.settled_start({k: v["pass"] for k, v in pairs.items()}) @@ -105,18 +125,29 @@ def main(argv): if "--reference" in argv: i = argv.index("--reference") ref_path, argv = argv[i + 1], argv[:i] + argv[i + 2:] + horizon = p.LAST_SCORED # default cap: bit-for-bit the original scoring + if "--horizon" in argv: + i = argv.index("--horizon") + horizon, argv = int(argv[i + 1]), argv[:i] + argv[i + 2:] paths = argv if ref_path is None or ref_path in argv else argv + [ref_path] cases = [load(path) for path in paths] - ledgers = [ledger(c) for c in cases] + for c in cases: # every scored cycle must be complete; never pad or shift + if c["cycles"] < horizon or len(c["force"]) < round((0.25 + horizon) * c["T"]): + raise SystemExit(f"{c['path']}: run does not contain complete cycles 1..{horizon}") + ledgers = [ledger(c, horizon) for c in cases] window = p.common_window([k for _, k in ledgers], - [{k: v["pass"] for k, v in pairs.items()} for pairs, _ in ledgers]) - report = dict(window=window, cases=[]) + [{k: v["pass"] for k, v in pairs.items()} for pairs, _ in ledgers], + last_scored=horizon) + report = dict(window=window, horizon=horizon, cases=[], + horizon_mode="default-cap" if horizon == p.LAST_SCORED else + "extended duration diagnostic (declared amendment, D40/60D only)") for case, (pairs, k_case) in zip(cases, ledgers): entry = {key: case[key] for key in ("path", "D", "T", "box", "n", "u_max", "nu", "tau", "A", "steps", "cycles", "F_ref")} entry.update(duration_T=case["steps"] / case["T"], K_case=k_case, pairs=pairs, coverage={k: coverage(case, k) for k in sorted(cycle_slices(case)) if k > p.DISCARD}, - no_donor=float(case["diag"][:, 13].sum()), eq_only=float(case["diag"][:, 12].sum())) + no_donor=float(case["diag"][:, 13].sum()), eq_only=float(case["diag"][:, 12].sum()), + hf_component=hf_component(case, window or tuple(range(horizon - 2, horizon + 1)))) if window: # residual drift inside the scored window (pairs K..K+2) entry["window_drift"] = {key: max(pairs[k][key] for k in window) for key in ("wave", "profile", "morison")} diff --git a/research/duetsch_re100_kc5/protocol.py b/research/duetsch_re100_kc5/protocol.py index be17314..1d65e49 100644 --- a/research/duetsch_re100_kc5/protocol.py +++ b/research/duetsch_re100_kc5/protocol.py @@ -114,14 +114,14 @@ def settled_start(pairs_pass): return None -def common_window(k_cases, pairs_pass_by_case): +def common_window(k_cases, pairs_pass_by_case, last_scored=LAST_SCORED): """Common scored cycles K..K+2 with K = max K_case, or None (FAIL-CLOSED) when a case never settled, K + 2 exceeds the cap, or a case fails any pair p_j for j = its own K_case .. K + 2 (continuous settling, including the wait for K).""" if any(k is None for k in k_cases): return None K = max(k_cases) - if K + 2 > LAST_SCORED: + if K + 2 > last_scored: return None for k_case, pairs in zip(k_cases, pairs_pass_by_case): if not all(pairs.get(j, False) for j in range(k_case, K + 3)): diff --git a/tests/unit/test_duetsch_protocol.py b/tests/unit/test_duetsch_protocol.py index 375411c..945633f 100644 --- a/tests/unit/test_duetsch_protocol.py +++ b/tests/unit/test_duetsch_protocol.py @@ -40,6 +40,10 @@ def test_settling_and_common_window(): assert p.common_window([7, 7], [passes, {**passes, 8: False}]) is None # must stay settled # a case settled early must keep passing while waiting for the common K assert p.common_window([7, 9], [{**passes, 8: False}, passes]) is None + long = {k: True for k in range(3, 32)} + assert p.common_window([10], [long]) is None # default cap 11 + assert p.common_window([29], [long], last_scored=31) == (29, 30, 31) + assert p.common_window([30], [long], last_scored=31) is None def test_solid_aware_bilinear():