From 70c83b7233888cb7287c8d1905563590c9e72ef4 Mon Sep 17 00:00:00 2001 From: TRON Date: Sun, 18 Oct 2026 21:25:43 -0700 Subject: [PATCH] test: add POST-G shear initialization counterfactual --- docs/validation/taylor_green.md | 28 ++++ tests/validation/taylor_green/run_case.py | 50 +++++-- tests/validation/taylor_green/shear_postg.py | 34 +++++ .../taylor_green/test_shear_wave.py | 129 +++++++++++++++++- 4 files changed, 225 insertions(+), 16 deletions(-) create mode 100644 tests/validation/taylor_green/shear_postg.py diff --git a/docs/validation/taylor_green.md b/docs/validation/taylor_green.md index 7703e99..61bd23e 100644 --- a/docs/validation/taylor_green.md +++ b/docs/validation/taylor_green.md @@ -480,3 +480,31 @@ can support a descriptive timing contrast, but cannot establish a crossing or reflection coefficient. Changed-level probes must not be treated as collocated samples. `rho-1` remains raw density deviation; do not substitute the 24-probe mean for the global pressure gauge. + +### POST-G initialization counterfactual (`postg_shear_n128`) + +Exactly three fresh Re100 cases are admitted: U128, U512, and projected R128 +depth 3 with the original nested rectangles. Each requires `--flow shear_wave +--diagnostic postg_shear_n128 --probe-trace`; snapshots and early-pressure mode +are excluded. Geometry, physical velocity/density, viscosities, collisions, +five main times and the 24 signed probes/windows are unchanged. Total work is +2,952,265,728 active updates, with 15 main rows, 774 trace rows and five builds. + +Only initialization changes. At level spacing h, k_level = h k and +A_i = 3 w_i U cos(theta) (c_i dot e) (c_i dot k_level). Initialize physical +FP64 **POST G = feq - (tau - 1) A**, then use the existing shifted FP32 codec. +This is the first-gradient, low-Mach Chapman-Enskog approximation, not an +exact discrete eigenstate. The PRE coefficient is -tau; the physical viscous +stress coefficient tau - 1/2 must not replace the raw POST-G coefficient. +The correction has zero mass/momentum before encoding; FP32 introduces normal +roundoff. U512 uses its native k and tau, matching nested L2 physical spacing. + +Receipts pin the new initializer module and identify `postg_first_gradient_v1`, +with `qualification_eligible=false`. Existing scalar screens remain visible +as descriptive outcomes, including failures; no new pass gate or assessor is +introduced. Historical assessors reject the new diagnostic and initialization +identities. Compare raw rows and signed traces with the retained equilibrium +receipts as a startup-state counterfactual. The host first-interval pressure +RMS reduction is not a prediction of long-time improvement, and individual +probe errors can increase. All equilibrium paths and historical verdicts stay +unchanged. diff --git a/tests/validation/taylor_green/run_case.py b/tests/validation/taylor_green/run_case.py index af36055..9d580c1 100644 --- a/tests/validation/taylor_green/run_case.py +++ b/tests/validation/taylor_green/run_case.py @@ -27,6 +27,8 @@ TWIN_SHEAR = "twin_shear_n128" # depth-2 twin: nested L0->L1 boxes without L2 SHEAR_TWIN_SHA256 = "1cbf37f18e89aea577581631c4565f51e68d22a3739b698816ef64acbf2b2bb4" WIDE_X_SHEAR = "wide_x_shear_n128" # one projected x-gap-width contrast SHEAR_WIDE_X_SHA256 = "13d16ba355117e2d074edccd9dc5d694c2fc13c673a0abdfca6401f6ed917733" +POSTG_SHEAR = "postg_shear_n128" # initialization-only three-case counterfactual +SHEAR_POSTG_SHA256 = "9c909fa013c092490ee21532e3125b2a78a4c5102c7a1bdc0e76986d3f0ef903" SHEAR_PROBES_SHA256 = "e2e6490d281cc275f448dc6d6ea7b9b4407cc74a51e9256b86c08c8626b89b76" # signed interface-probe mechanism diagnostic @@ -123,6 +125,13 @@ def verify_source(repo, revision, flow="taylor_green", probes=False): or hashlib.sha256(module.read_bytes()).hexdigest() != SHEAR_PROBES_SHA256): raise RuntimeError("Wrong shear_probes module checkout or content") result["shear_probes_sha256"] = SHEAR_PROBES_SHA256 + if flow == "shear_postg": + import shear_postg + module = harness / "shear_postg.py" + if (Path(shear_postg.__file__).resolve() != module + or hashlib.sha256(module.read_bytes()).hexdigest() != SHEAR_POSTG_SHA256): + raise RuntimeError("Wrong shear_postg module checkout or content") + result["shear_postg_sha256"] = SHEAR_POSTG_SHA256 if flow == "shear_wide_x": import shear_wide_x module = harness / "shear_wide_x.py" @@ -137,14 +146,14 @@ def verify_source(repo, revision, flow="taylor_green", probes=False): or hashlib.sha256(module.read_bytes()).hexdigest() != SHEAR_TWIN_SHA256): raise RuntimeError("Wrong shear_twin module checkout or content") result["shear_twin_sha256"] = SHEAR_TWIN_SHA256 - if flow in ("shear_nested", "shear_twin", "shear_wide_x"): # shared nested constants + if flow in ("shear_nested", "shear_twin", "shear_wide_x", "shear_postg"): import shear_nested module = harness / "shear_nested.py" if (Path(shear_nested.__file__).resolve() != module or hashlib.sha256(module.read_bytes()).hexdigest() != SHEAR_NESTED_SHA256): raise RuntimeError("Wrong shear_nested module checkout or content") result["shear_nested_sha256"] = SHEAR_NESTED_SHA256 - if flow in ("shear_wave", "shear_nested", "shear_twin", "shear_wide_x"): + if flow in ("shear_wave", "shear_nested", "shear_twin", "shear_wide_x", "shear_postg"): import shear_wave module = harness / "shear_wave.py" if (Path(shear_wave.__file__).resolve() != module @@ -218,15 +227,18 @@ def main(): p.add_argument("--output", type=Path, required=True) p.add_argument("--explosion", choices=["constant", "projected_linear"], default="constant") p.add_argument("--diagnostic", - choices=[SCALED_DIAGNOSTIC, SCALED_SHEAR, NESTED_SHEAR, TWIN_SHEAR, WIDE_X_SHEAR], + choices=[SCALED_DIAGNOSTIC, SCALED_SHEAR, NESTED_SHEAR, TWIN_SHEAR, + WIDE_X_SHEAR, POSTG_SHEAR], help="explicit four-case N256/fine512 same-placement diagnostic " "(scaled_n256: Taylor-Green; scaled_shear_n256: shear wave) or the " "fixed three-level/two-region shear preset (nested_shear_n128) " "or its depth-2 twin (twin_shear_n128), or the single projected " - "x-gap-width contrast (wide_x_shear_n128)") + "x-gap-width contrast (wide_x_shear_n128), or three-case POST-G " + "initialization counterfactual (postg_shear_n128)") p.add_argument("--probe-trace", action="store_true", help="signed interface-probe trace (nested U128/U512/projected d3, " - "twin projected d2, wide-x projected d3); mechanism diagnostic, no gates") + "twin projected d2, wide-x projected d3, POST-G three-case set); " + "mechanism diagnostic, no gates") p.add_argument("--final-snapshot", action="store_true", help="nested_shear_n128 projected depth-3 only: final macro snapshot") p.add_argument("--early-pressure", action="store_true", @@ -235,18 +247,21 @@ def main(): args = p.parse_args() nested = args.diagnostic == NESTED_SHEAR wide_x = args.diagnostic == WIDE_X_SHEAR - if wide_x and not args.probe_trace: - p.error("wide_x_shear_n128 requires --probe-trace") - if args.flow == "shear_wave" and (args.re != 100 or (args.depth == 3 and not (nested or wide_x)) + postg = args.diagnostic == POSTG_SHEAR + if (wide_x or postg) and not args.probe_trace: + p.error(f"{args.diagnostic} requires --probe-trace") + if args.flow == "shear_wave" and (args.re != 100 or (args.depth == 3 and not (nested or wide_x or postg)) or args.early_pressure): p.error("shear_wave requires --re 100, depth 1 or 2 and no --early-pressure") probe_cases = {(NESTED_SHEAR, 128, 1, "constant"), (NESTED_SHEAR, 512, 1, "constant"), (NESTED_SHEAR, 128, 3, "projected_linear"), (TWIN_SHEAR, 128, 2, "projected_linear"), - (WIDE_X_SHEAR, 128, 3, "projected_linear")} + (WIDE_X_SHEAR, 128, 3, "projected_linear"), + (POSTG_SHEAR, 128, 1, "constant"), (POSTG_SHEAR, 512, 1, "constant"), + (POSTG_SHEAR, 128, 3, "projected_linear")} if args.probe_trace and (args.final_snapshot or ( args.diagnostic, args.n, args.depth, args.explosion) not in probe_cases): - p.error("--probe-trace requires a fixed nested, twin or wide-x case, " + p.error("--probe-trace requires a fixed nested, twin, wide-x or POST-G case, " "without --final-snapshot") if args.final_snapshot and not (nested and args.depth == 3 and args.explosion == "projected_linear"): @@ -254,14 +269,18 @@ def main(): initial, metrics, limits, initialization = flow_functions(args.flow) if args.output.exists(): p.error("Output already exists; choose a new path to retain prior evidence") - if args.explosion != "constant" and args.depth != 2 and not (nested or wide_x): + if args.explosion != "constant" and args.depth != 2 and not (nested or wide_x or postg): p.error("projected_linear requires depth 2") twin = args.diagnostic == TWIN_SHEAR - if nested or twin or wide_x: + if nested or twin or wide_x or postg: import shear_nested import shear_twin import shear_wide_x preset = shear_wide_x if wide_x else shear_nested if nested else shear_twin + if postg: + import shear_postg + preset = shear_postg + initial, initialization = preset.initial, preset.INITIALIZATION if args.flow != "shear_wave": p.error(f"--diagnostic {args.diagnostic} requires --flow shear_wave") if (args.re, args.n, args.depth, args.explosion) not in preset.CASES: @@ -285,7 +304,7 @@ def main(): p.error("--run-gpu and exactly one explicit CUDA_VISIBLE_DEVICES are required") repo = args.repo.resolve() sys.path[:0] = [str(repo / "src"), str(repo)] - source_flow = ("shear_wide_x" if wide_x else "shear_nested" if nested + source_flow = ("shear_postg" if postg else "shear_wide_x" if wide_x else "shear_nested" if nested else "shear_twin" if twin else args.flow) source = verify_source(repo, args.expected_revision, source_flow, probes=args.probe_trace) runtime = runtime_source(repo) @@ -303,7 +322,7 @@ def main(): cfg = config(args.n, args.re) if args.depth == 1: refinements = [] - elif nested or twin or wide_x: + elif nested or twin or wide_x or postg: refinements = [list(level) for level in preset.REFINEMENTS] else: refinements = rectangles(args.n, args.depth) @@ -351,6 +370,9 @@ def main(): phase=shear_wave.PHI, eps=shear_wave.eps(args.n)) if args.diagnostic: report["diagnostic"] = args.diagnostic + if postg: + report.update(initialization_id=preset.INITIALIZATION_ID, qualification_eligible=False, + screening_scope="descriptive scalar screens only; no qualification assessor") args.output.parent.mkdir(parents=True, exist_ok=True) def save(): temporary = args.output.with_suffix(args.output.suffix + ".tmp") diff --git a/tests/validation/taylor_green/shear_postg.py b/tests/validation/taylor_green/shear_postg.py new file mode 100644 index 0000000..fd1af1a --- /dev/null +++ b/tests/validation/taylor_green/shear_postg.py @@ -0,0 +1,34 @@ +"""First-gradient POST-G initialization counterfactual, not a qualification preset. + +PRE F = feq - tau*A; POST G = feq - (tau - 1)*A. This low-Mach, +first-order Chapman-Enskog approximation is not an exact discrete eigenstate. +The physical viscous-stress factor tau - 1/2 is not the raw POST-G factor. +""" + +import numpy as np + +import shear_nested as nested +import shear_wave as sw +from oracle import CX, CY, W, U, centers + +REFINEMENTS = nested.REFINEMENTS +CASES = ((100, 128, 1, "constant"), (100, 512, 1, "constant"), + (100, 128, 3, "projected_linear")) +INITIALIZATION_ID = "postg_first_gradient_v1" +INITIALIZATION = ( + "physical FP64 POST G=feq-(tau-1)A at actual cell centers; " + "A=3*w*U*cos(theta)*(c.e)*(c.k_level); rho=1; first-gradient low-Mach CE" +) + + +def initial(shape, level, length, viscosity): + """Return physical FP64 G; production initialize owns shifted FP32 encoding.""" + x, y = centers(shape, level) + k = sw.wavevector(length) + theta = k[0] * x + k[1] * y + sw.PHI + k_level = k * 2.0**-level + tau = 0.5 + 3 * viscosity * 2**level + ce = CX * sw.E[0] + CY * sw.E[1] + ck = CX * k_level[0] + CY * k_level[1] + correction = (3 * W * ce * ck)[:, None, None] * U * np.cos(theta) + return sw.initial(shape, level, length, viscosity) - (tau - 1) * correction diff --git a/tests/validation/taylor_green/test_shear_wave.py b/tests/validation/taylor_green/test_shear_wave.py index f828450..05c5bb5 100644 --- a/tests/validation/taylor_green/test_shear_wave.py +++ b/tests/validation/taylor_green/test_shear_wave.py @@ -1062,7 +1062,7 @@ class ProbeTraceTests(unittest.TestCase): with self.assertRaises(SystemExit): run_case.main() - def run_mocked(self, n, depth, diagnostic, scheme, masks): + def run_mocked(self, n, depth, diagnostic, scheme, masks, initializer=sw.initial): from types import SimpleNamespace from contextlib import nullcontext arrays = self.exact_arrays(masks, n, 0.0) @@ -1081,7 +1081,9 @@ class ProbeTraceTests(unittest.TestCase): pass def initialize(self, *given): - pass + assert len(given) == depth + for level, (values, mask) in enumerate(zip(given, masks)): + np.testing.assert_array_equal(values, initializer(mask.shape, level, n, U * n / sw.RE)) def step(self, count): self.step_count += count @@ -1289,3 +1291,126 @@ class WideXDiagnosticTests(unittest.TestCase): refinements=self.wide.REFINEMENTS) self.assertEqual(helper.mod.assess_nested(records, "r")["status"], "INVALID_SET") self.assertNotIn(run_case.WIDE_X_SHEAR, helper.mod.PRESETS) + + +class PostGDiagnosticTests(unittest.TestCase): + """Three initialization-only cases, excluded from historical qualification.""" + + def setUp(self): + import shear_postg + self.postg = shear_postg + + def test_moments_phase_and_physical_scaling(self): + from CelerisLab.lbm.population_storage import _encode_physical, _decode_stored + velocities = np.array([CX, CY], dtype=float).T + for level in range(3): + shape, nu = (31, 47), .0384 + h, tau = 2.0**-level, .5 + 3 * nu * 2**level + eq = sw.initial(shape, level, 128, nu) + g = self.postg.initial(shape, level, 128, nu) + x, y = centers(shape, level) + k = sw.wavevector(128) + strain = (.5 * U * np.cos(k[0] * x + k[1] * y + sw.PHI)[None, None] + * (np.outer(sw.E, h * k) + np.outer(h * k, sw.E))[:, :, None, None]) + tensor = np.einsum("qi,qj->qij", velocities, velocities) - np.eye(2)[None] / 3 + a = 3 * oracle.W[:, None, None] * np.einsum("qij,ijyx->qyx", tensor, strain) + np.testing.assert_allclose(g, eq - (tau - 1) * a, rtol=0, atol=6e-17) + for weights in (np.ones(9), CX, CY): + np.testing.assert_allclose(np.einsum("q,qyx->yx", weights, g - eq), + 0, rtol=0, atol=2e-16) + stress = np.einsum("qi,qj,qyx->ijyx", velocities, velocities, g - eq) + np.testing.assert_allclose(stress, (1 - tau) * 2 / 3 * strain, rtol=0, atol=6e-17) + pre = eq - tau * a + rho = pre.sum(axis=0) + pre_eq = sw.equilibrium(rho, np.einsum("q,qyx->yx", CX, pre) / rho, + np.einsum("q,qyx->yx", CY, pre) / rho) + np.testing.assert_allclose((1 - 1 / tau) * pre + pre_eq / tau, g, rtol=0, atol=4e-16) + self.assertAlmostEqual((tau - .5) * h, .1152, places=16) + self.assertGreater(1 - tau, 0) + encoded = _decode_stored(_encode_physical(g, True), True) + self.assertGreater(encoded.min(), 0) + for weights in (np.ones(9), CX, CY): + np.testing.assert_allclose(np.einsum("q,qyx->yx", weights, encoded - eq), + 0, rtol=0, atol=3e-9) + # The U512 control is the same physical spacing/tau/phase as nested L2. + np.testing.assert_array_equal(self.postg.initial((512, 512), 0, 512, .1536), + self.postg.initial((512, 512), 2, 128, .0384)) + + def test_actual_initializers_and_three_complete_schedules(self): + helper = ProbeTraceTests() + helper.setUp() + work = 0 + for _, n, depth, scheme in self.postg.CASES: + masks = helper.nested.leaf_masks() if depth == 3 else [np.ones((n, n), bool)] + rec, counts, files = helper.run_mocked(n, depth, run_case.POSTG_SHEAR, scheme, masks, + initializer=self.postg.initial) + self.assertEqual(rec["initialization"], self.postg.INITIALIZATION) + self.assertEqual(rec["initialization_id"], "postg_first_gradient_v1") + self.assertFalse(rec["qualification_eligible"]) + self.assertTrue(rec["completed"]) + self.assertEqual(rec["limits"], sw.LIMITS) + self.assertEqual([r["step"] for r in rec["rows"]], run_case.sample_steps(n)[0]) + self.assertEqual(rec["probe_trace"]["steps"], list(helper.probes.STEPS)) + self.assertEqual(rec["probe_trace"]["probes"], helper.probes.select(masks, n)) + self.assertEqual(counts[-1], 16 * n) + self.assertIn(1920 * (n // 128), counts) + self.assertEqual(files, ["probe.json"]) + self.assertNotIn("final_snapshot", rec) + work += rec["workload"]["active_collision_updates"] + self.assertEqual(work, 2952265728) + + def test_fixed_admission_and_source_pin(self): + from tests.validation.multigrid import check_projected_pairs as gate + repo = HERE.parents[2] + self.assertEqual(run_case.SHEAR_POSTG_SHA256, + hashlib.sha256((HERE / "shear_postg.py").read_bytes()).hexdigest()) + with tempfile.TemporaryDirectory() as tmp: + base = ["run_case.py", "--repo", str(repo), "--expected-revision", "f" * 40, + "--output", str(Path(tmp) / "case.json"), "--dry-run", "--re", "100", + "--flow", "shear_wave", "--diagnostic", run_case.POSTG_SHEAR, "--probe-trace"] + for n, depth, scheme, valid in ((128, 1, "constant", True), (512, 1, "constant", True), + (128, 3, "projected_linear", True), + (128, 3, "constant", False), (256, 1, "constant", False), + (128, 2, "projected_linear", False)): + argv = base + ["--n", str(n), "--depth", str(depth), "--explosion", scheme] + extras = ([], ["--final-snapshot"], ["--early-pressure"], ["--re", "1000"], + ["--flow", "taylor_green"]) if valid else ([],) + for extra in extras: + with self.subTest(argv=argv + extra), patch.object(sys, "argv", argv + extra), \ + patch.object(run_case, "verify_source", return_value={}) as verify, \ + patch.object(run_case, "runtime_source", return_value={}), \ + redirect_stdout(io.StringIO()), patch("sys.stderr", io.StringIO()): + if valid and not extra: + self.assertEqual(run_case.main(), 0) + verify.assert_called_once_with(repo, "f" * 40, "shear_postg", probes=True) + else: + with self.assertRaises(SystemExit): + run_case.main() + verify.assert_not_called() + with patch.object(sys, "argv", base[:-1] + ["--n", "128", "--depth", "1"]), \ + patch("sys.stderr", io.StringIO()): + with self.assertRaises(SystemExit): + run_case.main() + with patch.object(gate, "preflight", return_value=run_case.NUMERICAL_REVISION), \ + patch.object(run_case, "provenance", return_value=dict(source_tree_sha256=run_case.SOURCE_SHA256)): + self.assertEqual(run_case.verify_source(repo, "r", "shear_postg")["shear_postg_sha256"], + run_case.SHEAR_POSTG_SHA256) + for module, attr, value in ((self.postg, "__file__", "/wrong/shear_postg.py"), + (run_case, "SHEAR_POSTG_SHA256", "wrong")): + with patch.object(module, attr, value), self.assertRaisesRegex(RuntimeError, "shear_postg"): + run_case.verify_source(repo, "r", "shear_postg") + + def test_old_assessors_reject_identity_or_initializer(self): + from copy import deepcopy + helper = NestedPresetTests() + helper.setUp() + for base, assessor in ((helper.nested_set(), helper.mod.assess_nested), + ([receipt(128, 1), receipt(256, 1), receipt(128, 2), + receipt(128, 2, "projected_linear")], helper.mod.assess)): + self.assertTrue(assessor(base, "r")["passed"]) + for change in (dict(diagnostic=run_case.POSTG_SHEAR), + dict(initialization=self.postg.INITIALIZATION)): + records = deepcopy(base) + records[-1].update(change) + self.assertEqual(assessor(records, "r")["status"], "INVALID_SET") + self.assertNotIn(run_case.POSTG_SHEAR, helper.mod.PRESETS)