test: DFG Re20 steady rule per review (1 % band, range and slope, 0.5 s samples)

Samples every 0.5 s at synchronized root steps; trailing 5 s window (11
samples); range and |slope*window| each at most 1 % of the official band
width. States the fluid-on-cylinder force sign, reports cd_positive,
defaults to the refined case and uses the review's label.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>

docs/validation/dfg_re20.md

@@ -6,7 +6,7 @@ inflow, initial state, layout and observation are those of the pinned Re100 runn
 (`run_refined_dfg.py`, `re100_helpers.py`, both unchanged); only the lattice
 viscosity is set for Re 20 (root tau 0.56, fine 0.62, uniform 0.62, lattice speed 0.01).
 
-This is a comparison with the official intervals, not a full benchmark
+Label: "compared against the official reference bands", not a full benchmark
 certification: the outlet is pressure-NEQ (`neq_extrap`), not do-nothing, and the
 inlet is local Zou-He.
 
@@ -22,10 +22,13 @@ inlet is local Zou-He.
 ## Bounds and steady rule (frozen before any run)
 
 - End time t* = 60 (30 s physical): 240,000 refined root steps, 480,000 uniform steps.
-- One sample every t* = 0.5; trailing window t* = 10 (5 s, 20 samples).
+- One sample every t* = 1 (0.5 s) at the end of a synchronized root step; trailing
-- Steady: in the window, max - min of each of Cd, Cl, p_diff is at most one tenth of
+  window t* = 10 (last 5 s, 11 samples).
-  its official interval width (0.002, 0.00006, 0.00004). Absolute ranges, so the
+- Steady: in the window, both max - min and |least-squares slope x window| of each
-  small Cl is never used as a divisor.
+  of Cd, Cl, p_diff are at most 1 % of its official interval width (2e-4, 6e-6,
+  4e-6). Absolute limits, so the small Cl is never a divisor. Never relaxed after a run.
+- Force sign: fluid on cylinder, x along the inflow, y up from the bottom wall; Cd > 0
+  and the official Cl > 0 (`cd_positive` is reported).
 - Not steady at t* = 60: status `NOT_STEADY`, no interval verdict. Health limits
   (Mach <= 0.2, rho in [0.5, 1.5]) stop the run with `UNHEALTHY`.
 - Steady: `WITHIN_INTERVALS` or `OUTSIDE_INTERVALS`, with the per-quantity verdict.
@@ -34,10 +37,10 @@ inlet is local Zou-He.
 
 ```
 python tests/validation/run_dfg_re20.py --output <new-dir> --plan
-CUDA_VISIBLE_DEVICES=<n> python tests/validation/run_dfg_re20.py --case both \
+CUDA_VISIBLE_DEVICES=<n> python tests/validation/run_dfg_re20.py --case refined \
     --collision SRT --output <new-dir> --max-wall-seconds <cap>
 ```
 
-`--plan` writes `plan.json` only. A GPU run writes `plan.json` and `result.json`
+Run `--case uniform` only if the refined result is outside a band. `--plan` writes `plan.json` only. A GPU run writes `plan.json` and `result.json`
 (per case: lattice parameters, every sample, window ranges, status, verdict).
 Host contract tests: `tests/unit/test_dfg_re20.py`.

tests/unit/test_dfg_re20.py

@@ -21,8 +21,8 @@ class Re20Contract(unittest.TestCase):
         root, uniform = v.lattice_parameters(40), v.lattice_parameters(80)
         self.assertAlmostEqual(root["tau"], 0.56, places=12)
         self.assertAlmostEqual(uniform["tau"], 0.62, places=12)
-        self.assertEqual((root["steps"], root["interval"]), (240000, 2000))
+        self.assertEqual((root["steps"], root["interval"]), (240000, 4000))
-        self.assertEqual((uniform["steps"], uniform["interval"]), (480000, 4000))
+        self.assertEqual((uniform["steps"], uniform["interval"]), (480000, 8000))
         self.assertAlmostEqual(v.plan()["refined_fine_tau"], 0.62, places=12)
 
     def test_config_only_changes_viscosity(self):
@@ -47,17 +47,22 @@ class Re20Contract(unittest.TestCase):
     def test_steady_and_verdict(self):
         inside = dict(cd=5.58, cl=0.0106, p_diff=0.1174)
         out = v.summarize(rows(lambda i: inside))
-        self.assertEqual((out["status"], out["steady"]), ("WITHIN_INTERVALS", True))
+        self.assertEqual((out["status"], out["steady"], out["window_samples"]),
+                         ("WITHIN_INTERVALS", True, 11))
+        self.assertTrue(out["cd_positive"])
         self.assertTrue(all(out["verdict"].values()))
         low = v.summarize(rows(lambda i: dict(inside, cl=0.0100)))
         self.assertEqual(low["status"], "OUTSIDE_INTERVALS")
         self.assertFalse(low["verdict"]["cl"])
 
     def test_drift_gives_no_verdict(self):
-        # Cl drifts by more than a tenth of its interval width inside the window.
+        # Cl drifts by more than 1 % of its interval width (6e-6) inside the window.
-        out = v.summarize(rows(lambda i: dict(cd=5.58, cl=0.0106 + 1e-5 * i, p_diff=0.1174)))
+        out = v.summarize(rows(lambda i: dict(cd=5.58, cl=0.0106 + 1e-6 * i, p_diff=0.1174)))
         self.assertEqual((out["status"], out["verdict"]), ("NOT_STEADY", None))
         self.assertGreater(out["ranges"]["cl"], out["steady_limits"]["cl"])
+        self.assertGreater(out["drifts"]["cl"], out["steady_limits"]["cl"])
+        small = v.summarize(rows(lambda i: dict(cd=5.58, cl=0.0106 + 2e-6 * (i % 2), p_diff=0.1174)))
+        self.assertTrue(small["steady"])  # 2e-6 alternation, no trend: within both limits
 
     def test_incomplete_run_gives_no_verdict(self):
         out = v.summarize(rows(lambda i: dict(cd=5.58, cl=0.0106, p_diff=0.1174), tstar_end=30))

tests/validation/run_dfg_re20.py

@@ -11,16 +11,24 @@ Cd = 2 Fx / (rho U^2 D) and Cl = 2 Fy / (rho U^2 D) are scale-free, and
 ``pressure_probes`` returns dp / (rho U^2). The official p_diff = p(0.15, 0.2) -
 p(0.25, 0.2) uses U = 0.2 m/s, rho = 1, so p_diff = dp * 0.04. Both are reported;
 only p_diff is compared with the official interval. Forces use the existing sign
-convention (positive Cd along the inflow, positive Cl towards +y).
+convention: F is the force of the fluid on the cylinder, i.e. the refined raw
+exchange impulse (window ledger, ``total.impulse``) or the uniform ``read_force``,
+with x along the inflow and y up from the bottom wall. Cd must be positive, and the
+official Cl is positive (centre 0.20 below the centreline 0.205), so a flipped sign
+shows up as Cl outside its interval.
 
-Sampling. One sample every SAMPLE_TSTAR (t* = U t / D): Cd and Cl are the mean
+Sampling. One sample every SAMPLE_TSTAR (t* = U t / D; 0.5 s), taken at the end of a
-over that interval (window exchange ledger refined, per-batch forces uniform);
+completed root step with all levels synchronized; the step index is recorded. Cd
-p_diff is the end-of-interval value at the exact surface points.
+and Cl are the mean over the sample interval (smooths step-parity/acoustic
+oscillation; equal to the instantaneous value at steady state); p_diff is the
+end-of-interval value at the exact surface points.
 
-Steady rule (frozen before any run). Over the trailing WINDOW_TSTAR, the range
+Steady rule (frozen before any run). Over the trailing WINDOW_TSTAR (last 5 s, 11
-(max - min) of each of Cd, Cl and p_diff must not exceed one tenth of the width of
+samples), both the range (max - min) and |least-squares slope * window| of each of
-its official interval. Absolute ranges avoid dividing by the small Cl. If the rule
+Cd, Cl and p_diff must not exceed 1 % of the width of its official interval. These
-fails at TSTAR_MAX the run reports "not steady" and gives no interval verdict.
+are absolute limits, so the small Cl is never a divisor, and the slope term catches
+slow monotone drift. If the rule fails at TSTAR_MAX the run reports "not steady" and
+gives no interval verdict; the limits are never relaxed after a run.
 """
 
 import argparse
@@ -41,11 +49,11 @@ import run_refined_dfg as rr  # noqa: E402
 RE = 20.0
 U_PHYS = 0.2                   # m/s mean inflow; rho = 1
 TSTAR_MAX = 60.0               # 30 s physical
-SAMPLE_TSTAR = 0.5
+SAMPLE_TSTAR = 1.0             # 0.5 s
 WINDOW_TSTAR = 10.0            # trailing 5 s
 # Schaefer & Turek (1996), 2D-1 lower/upper bounds.
 BANDS = dict(cd=(5.57, 5.59), cl=(0.0104, 0.0110), p_diff=(0.1172, 0.1176))
-STEADY_RANGE = {k: (hi - lo) / 10 for k, (lo, hi) in BANDS.items()}
+STEADY_LIMIT = {k: (hi - lo) / 100 for k, (lo, hi) in BANDS.items()}
 DEVIATIONS = ("pressure-NEQ outlet (neq_extrap), not do-nothing", "local Zou-He inlet")
 
 _re100_config = rr.config
@@ -78,11 +86,14 @@ def summarize(rows):
     if not rows or rows[-1]["tstar"] < TSTAR_MAX - 1e-9:
         return dict(status="INCOMPLETE", steady=False, verdict=None)
     window = [r for r in rows if r["tstar"] >= TSTAR_MAX - WINDOW_TSTAR - 1e-9]
+    t = np.array([r["tstar"] for r in window])
     ranges = {k: max(r[k] for r in window) - min(r[k] for r in window) for k in BANDS}
-    steady = all(ranges[k] <= STEADY_RANGE[k] for k in BANDS)
+    drifts = {k: abs(np.polyfit(t, [r[k] for r in window], 1)[0]) * WINDOW_TSTAR for k in BANDS}
+    steady = all(ranges[k] <= STEADY_LIMIT[k] and drifts[k] <= STEADY_LIMIT[k] for k in BANDS)
     final = {k: rows[-1][k] for k in BANDS}
-    result = dict(window_samples=len(window), ranges=ranges, steady_limits=STEADY_RANGE,
+    result = dict(window_samples=len(window), ranges=ranges, drifts=drifts,
-                  steady=steady, final=final, bands=BANDS, deviations=DEVIATIONS)
+                  steady_limits=STEADY_LIMIT, steady=steady, final=final, bands=BANDS,
+                  deviations=DEVIATIONS, cd_positive=final["cd"] > 0)
     if not steady:
         return dict(result, status="NOT_STEADY", verdict=None)
     inside = {k: BANDS[k][0] <= final[k] <= BANDS[k][1] for k in BANDS}
@@ -171,16 +182,16 @@ def run_case(a, h, case):
 
 def plan(collision="SRT"):
     return dict(re=RE, collision=collision, tstar_max=TSTAR_MAX, sample_tstar=SAMPLE_TSTAR,
-                window_tstar=WINDOW_TSTAR, bands=BANDS, steady_limits=STEADY_RANGE,
+                window_tstar=WINDOW_TSTAR, bands=BANDS, steady_limits=STEADY_LIMIT,
                 deviations=DEVIATIONS, pressure_scale=U_PHYS**2,
                 refined_root=lattice_parameters(40), refined_fine_tau=0.5 + 3 * 2 * rr.SPEED * 40 / RE,
                 uniform=lattice_parameters(80),
-                label="compared with official reference intervals; not a full benchmark certification")
+                label="compared against the official reference bands (non-canonical neq_extrap outlet, zou_he_local inlet); not a full benchmark certification")
 
 
 def main():
     parser = argparse.ArgumentParser(description=__doc__.splitlines()[0])
-    parser.add_argument("--case", choices=["both", "refined", "uniform"], default="both")
+    parser.add_argument("--case", choices=["refined", "uniform", "both"], default="refined")
     parser.add_argument("--collision", choices=["SRT", "MRT"], default="SRT")
     parser.add_argument("--source", type=Path, default=rr.DEFAULT_SOURCE)
     parser.add_argument("--output", type=Path, required=True)