From 413f7afade511cfcfd5758182e548da3b0ba9c06 Mon Sep 17 00:00:00 2001 From: TRON Date: Sun, 18 Oct 2026 00:15:52 -0700 Subject: [PATCH] test: integrate matched refined DFG pilot preparation --- README.md | 1 + docs/validation/refined_dfg.md | 70 ++ tests/unit/test_refined_dfg.py | 293 +++++++++ .../data/refined_dfg/source_manifest.json | 110 ++++ tests/validation/run_refined_dfg.py | 603 ++++++++++++++++++ 5 files changed, 1077 insertions(+) create mode 100644 docs/validation/refined_dfg.md create mode 100644 tests/unit/test_refined_dfg.py create mode 100644 tests/validation/data/refined_dfg/source_manifest.json create mode 100644 tests/validation/run_refined_dfg.py diff --git a/README.md b/README.md index be95081..191cb5a 100644 --- a/README.md +++ b/README.md @@ -665,6 +665,7 @@ streak.render("streakline.png") ## Validation campaigns +- [Matched refined DFG](docs/validation/refined_dfg.md): frozen SRT uniform/refined comparison, requiring channel GPU gates before execution. - [DFG 2D-2 Re100](docs/validation/dfg_re100.md): frozen periodicity and official-reference comparisons, with a declared noncanonical outlet. - [Uniform LES shear decay](docs/validation/uniform_les_flow.md): independent continuum-oracle flow diagnostic, baseline and isolated corrected source variants. - [Smagorinsky helper comparison](docs/validation/les_helpers.md): pinned baseline versus an output-local two-coefficient candidate, without production edits. diff --git a/docs/validation/refined_dfg.md b/docs/validation/refined_dfg.md new file mode 100644 index 0000000..375c628 --- /dev/null +++ b/docs/validation/refined_dfg.md @@ -0,0 +1,70 @@ +# Matched SRT uniform / refined DFG Re100 preparation + +Status: CPU-preflight checked; **not GPU executed or approved to execute**. No physical acceptance or speedup claim. No production source changes. Root owns GPU admission and integration. + +## Frozen setup + +Canonical source snapshot: integrated commit `cd3a657f9d3833e812c3cca173525f165b48710a`. Production origin remains candidate `47b05bc8abfc957fe876b0c63bdb386e3c6cbaa5` (on `5338864`, on `69e9b78`); production files are byte-identical. The integrated Re100 helper header/PIN records that origin explicitly. Both commit identities are retained in the manifest. `tests/validation/data/refined_dfg/source_manifest.json` hashes all production source files plus the existing DFG geometry/pressure helper, periodicity helper and official reference cycle. An integrated checkout may have a different commit but must have identical pinned file contents. No MRT external-reference run is reused as the matched uniform control. + +Both cases: SRT, shifted FP32, double buffer, no LES, no applied force, stationary cylinder, local Zou–He parabolic inlet, halfway y walls, identical pressure-NEQ outlet with **all nine nonequilibrium components damped by 0.5**, backflow clamp disabled. This outlet is not the exact canonical DFG do-nothing condition. Raw cylinder reaction observation is a diagnostic, not applied forcing. + +Physical units use D=0.1 m, mean U=1 m/s, ν=0.001 m²/s, density reference=1; domain 2.2×0.41 m, cylinder center (0.2,0.2) m. Uniform D80: 1761×330, dx=0.00125 m, dt=0.0000125 s, ν_lattice=.008. Refined root D40: 881×166, dx=.0025 m, dt=.000025 s, ν_lattice=.004; fine dx/dt/ν match uniform D80. All lattice mean speeds are .01. + +Fixed half-open patch **(40,24,320,140)** in root edge coordinates maps to x=[.09875,.79875] m, y=[.0575,.3475] m. The body and about six diameters of downstream wake are fine; wall layers and far wake remain coarse. Bounds are frozen before outcomes, not tuned afterward. Geometry-only leaf-count and nominal collision-work ratios to uniform are ~.419 and ~.321. These omit real launch/transfer/halo/measurement costs and are not measured speedups. The narrow root block height is legal but potentially expensive; do not infer speed from leaf count. + +Exact root dimensions are tiled by (by,bx)=(2,881); no padding. Fine allocated shape (240,568), including four-cell halos, uses (8,8). `host_geometry` calls production channel transfer-overlap checks and complete fine cut-link/donor/transfer guards before any CUDA import. + +### Half-cell conversion + +Uniform integer nodes have physical coordinates (i*dx, (j-.5)*dx). Root edge coordinates therefore map to ((X-.5)*dx_root, (Y-1)*dx_root). Body edge center=(80.5,81), radius=20; wall edges Y=1,165; source planes X=.5,880.5. Fine local body nodes are `(body_edge - fine_origin)/.5 - .5`. + +Body and physical wall planes match exactly. The axes have different phases. With fine origin O=(38,22): + +- x_f = [O_x+(j+.5)/2-.5] dx40 = (j+75.5) dx80, versus x_uniform = i dx80. Fine x samples have a half-fine-cell (quarter-root-cell) offset. +- y_f = [O_y+(k+.5)/2-1] dx40 = (k+42.5) dx80, versus y_uniform = (l-.5) dx80. The y samples align exactly for l=k+43. + +This difference follows from the inlet's ghost-node plane versus the half-way y-wall planes. The host phase test checks both mappings. Compare interpolated observables, not pointwise arrays or bitwise equality; x-phase error, the different fine collision arithmetic, wall resolution, and coarse/fine transfer effects all contribute to discrepancies. No field-difference acceptance is claimed by this runner. Initial fluid state samples the same analytic parabolic equilibrium as uniform production initialization; root boundary ghosts start at rest. It is not asserted to be bitwise identical across different grids/arithmetic. + +## Observations and gates + +- Every AMR base step is individually followed by `read_raw_exchange(..., window='last_base_step')`: two fine substeps, root-lattice per-unit-depth impulse integrated over [n,n+1]. `step(n)` alone retains only its last base interval, so it cannot support averaging the whole batch. +- Uniform `Simulation.run(batch)` / `read_force` gives the mean over that batch; multiply by batch before accumulating. Cd/Cl divide the interval mean reaction by 0.5 U_lattice² D_lattice in each case. Force unit is ρ dx³/dt² per depth, so fine-uniform lattice force is twice root force for the same physical load. +- Force timestamps are observation-interval midpoints. Pressure uses physical decoded postcollision density at the interval end, extrapolated from fluid-only points one/two **fine** cells outside the cylinder using existing `pressure_probes`. Pressure difference is normalized by ρ_ref U_physical². Existing periodic helpers align pressure onto force cycles. +- Public AMR downloads exclude nonleaves. Fixed geometry masks, not state-dependent finiteness masks, check all fluid populations and sum ρ dx_level²; covered root cells, body, boundaries and fine halos are excluded. Corrupt NaN fluid cannot disappear from the mass check. Mass is per-unit-depth ρ_ref·m². Uniform uses its physical fluid flags and dx². Different sampling of the circular boundary means initial discrete fluid volume can differ slightly; compare temporal stability, not identical absolute mass. +- Flux samples root/uniform x-index 2 and -3, summed with dx and physical velocity conversion. Those planes differ by grid-scale offsets; cycle-average inlet/outlet balance is the gate, not direct pointwise flux equality. +- Dense measurement width Δt*=.02 in the pilot and for t*≥240, otherwise .5. Pilot ends t*=5; full ends t*=300. Same existing five-cycle periodicity, mass, flux, and reference gates are imported without modification. Both official-reference diagnostics and refined-versus-uniform cycle comparisons are written. `matched_case_pass` refers only to matched SRT agreement, not certification of canonical DFG or the broader suite. +- Safety stop: nonfinite state/reaction, Mach>.2, or density outside [.5,1.5]. GPU durations include required exchange reads and synchronization in stepping; full wall time also includes host measurements. Do not report pure-kernel throughput from these numbers. + +## Commands + +From the repository root, CPU only: + +```sh +PYTHONPATH=src python tests/unit/test_refined_dfg.py +python tests/validation/run_refined_dfg.py --mode plan --output /tmp/refined-dfg-plan-new +``` + +Fourteen CPU tests are retained in tests/unit/test_refined_dfg.py. Plan JSON is configuration output only; no preparation result files are committed. The initial attempt to use pytest was unavailable in this interpreter; tests use standard-library unittest and NumPy instead. + +Only after root channel GPU gates pass and root explicitly admits an idle physical GPU, root may run: + +```sh +CUDA_VISIBLE_DEVICES= timeout --signal=TERM --kill-after=30s s python tests/validation/run_refined_dfg.py --mode pilot --gpu-gates-approved --max-wall-seconds --output +CUDA_VISIBLE_DEVICES= timeout --signal=TERM --kill-after=30s s python tests/validation/run_refined_dfg.py --mode full --gpu-gates-approved --max-wall-seconds --pilot-results --output +``` + +`--gpu-gates-approved` is a caller attestation, not proof that GPU tests ran. Full mode additionally requires both matched healthy pilots through t*=5 with identical source hashes, bounds and runner hash. No GPU/Mac/remote operation has been performed for this preparation. Pilot/full outputs must be new/empty. Uniform runs first; an unhealthy case stops the campaign. A full run that fails periodic/matched gates exits nonzero. + +## Bounded execution and receipt validation + +No pilot throughput is known yet; the caller must set a finite positive `--max-wall-seconds` rather than treating an invented estimate as calibration. After valid pilots, the full-run calibrated budget is 3 × sum_over_cases[startup_seconds + 60 × (pilot_wall_seconds − startup_seconds)]. The full run uses the smaller of that estimate and the caller's explicit ceiling. This uses total measured pilot evolution/observation cost and the full/pilot physical-time ratio (300/5); 3× is a declared safety margin, not a runtime guarantee. The plan records both budgets. + +The campaign deadline begins before the first case initialization. Cooperative checks surround initialization, every AMR root step or uniform ≤100-step batch, and measurements. On observed expiry, the runner stops, writes a nonpassing `wall_time_limit` receipt, and exits 3. A blocking CUDA call cannot be interrupted by these checks: the external `timeout` watchdog is required, with a hard budget chosen slightly above the soft budget to permit cleanup. Hard termination may leave partial CSV data and no final receipt; that never qualifies as a completed pilot. + +Each case emits `READY` with strict finite JSON and saves `ready.json` only after its first healthy complete observation, after initialization, actual stepping, force readback, macro/pressure/mass/flux paths and CSV flush all succeeded. No early initialization-only READY is emitted. + +Full-run pilot validation rejects nonfinite JSON constants and numeric overflow, wrong case/mode/status, nonboolean health, missing/duplicate/out-of-order rows, wrong step counts or observation times, missing numeric observations, health-limit violations, absent/mismatched READY receipts and invalid timing calibration. It requires exactly 250 observations spanning the full dense t*=5 pilot, source/runner hashes and frozen bounds. Tests include the reviewer's exact NaN-endpoint regression and malformed receipt variants. Synthetic test fixtures exist only in temporary directories and are not simulation results. + +## Integrated source import guard + +The runner's source defaults to this repository root. Before geometry or GPU work, it verifies the pinned source/helper/reference contents and requires the imported CelerisLab package to be exactly `/src/CelerisLab/__init__.py`. An installed DynamisLab copy or another cached checkout is rejected. Integration changes paths/readability and this explicit import assertion; geometry, numerical definitions, gates, and receipt validation remain frozen. Any new path-formatted runner has its own hash, so old preparation pilot receipts cannot authorize a full run of this integrated runner. diff --git a/tests/unit/test_refined_dfg.py b/tests/unit/test_refined_dfg.py new file mode 100644 index 0000000..a5c6f4e --- /dev/null +++ b/tests/unit/test_refined_dfg.py @@ -0,0 +1,293 @@ +"""CPU-only contract tests. These are not CUDA or DFG acceptance evidence.""" + +import json +import sys +from pathlib import Path +import numpy as np +import unittest + +sys.path.insert(0, str(Path(__file__).resolve().parents[1] / "validation")) +import run_refined_dfg as r + +h, periodic = r.helpers(r.DEFAULT_SOURCE) + + +def test_source_pin(): + assert r.verify_source(r.DEFAULT_SOURCE)["commit"].startswith("cd3a657") + + +def test_matched_units_and_planes(): + p = r.plan(h) + c = p["root"] + f = p["uniform"] + assert c["nx"] == 881 and c["ny"] == 166 + assert f["nx"] == 1761 and f["ny"] == 330 + assert c["dx"] / 2 == f["dx"] and c["dt"] / 2 == f["dt"] + assert np.isclose(c["viscosity"] * c["dx"] ** 2 / c["dt"], 0.001) + assert np.isclose(f["viscosity"] * f["dx"] ** 2 / f["dt"], 0.001) + assert np.isclose((p["body_root_edge"][0] - 0.5) * c["dx"], 0.2) + assert np.isclose((p["body_root_edge"][1] - 1) * c["dx"], 0.2) + assert np.isclose((165 - 1) * c["dx"], 0.41) + assert np.isclose((880.5 - 0.5) * c["dx"], 2.2) + assert np.isclose(2 * p["radius_root"] * c["dx"], 0.1) + + +def test_same_numerics(): + c, _, _ = r.config(h, 40) + f, _, _ = r.config(h, 80) + assert c["method"] == f["method"] + m = c["method"] + assert m["collision"] == "SRT" and m["streaming"] == "double_buffer" + assert m["store_precision"] == "FP32" and m["ddf_shifting"] + assert not m["les"]["enabled"] + assert m["inlet"]["scheme"] == "zou_he_local" + assert m["outlet"]["mode"] == "neq_extrap" and m["outlet"]["srt_neq_damp"] == 0.5 + assert not m["outlet"]["backflow_clamp"] + + +def test_actual_geometry_guards(): + cfg, pair, body, host = r.host_geometry(r.DEFAULT_SOURCE, h) + assert pair.fine_shape == (240, 568) + assert np.allclose( + host.objects[0].geometry.center, + [ + (80.5 - pair.fine_origin[0]) / 0.5 - 0.5, + (81 - pair.fine_origin[1]) / 0.5 - 0.5, + ], + ) + assert cfg.nx % r.BASE_BLOCK[1] == 0 and cfg.ny % r.BASE_BLOCK[0] == 0 + assert len(host.boundary_cells) > 0 + + +def test_physical_leaf_mass_excludes_covered_and_halo(): + c = np.full((9, 4, 4), 1 / 9) + f = np.full((9, 8, 8), 1 / 9) + cm = np.ones((4, 4), bool) + cm[1:3, 1:3] = False + fm = np.zeros((8, 8), bool) + fm[2:6, 2:6] = True + c[:, ~cm] = np.nan + f[:, ~fm] = np.nan + assert np.isclose(r.leaf_mass([c, f], 2, [cm, fm]), 64) + f[:, 3, 3] = np.nan + with unittest.TestCase().assertRaisesRegex(ValueError, "Nonfinite"): + r.leaf_mass([c, f], 2, [cm, fm]) + + +def test_force_impulse_interval_scaling(): + from CelerisLab.lbm.raw_exchange import summarize_exchange + + # Two fine substeps each j=(2,4), h=.5 => root impulse=(1,2). + result = summarize_exchange( + [[1.0, 2.0, 0.0, 0.0]], + [0], + [False], + body_id=0, + level=1, + h=0.5, + time_end=8, + window="last_base_step", + reference_center=(80.5, 81.0), + ) + assert result["sample_count"] == 2 and result["time_start"] == 7 + np.testing.assert_allclose(result["total"]["force"], [1, 2]) + assert result["units"] == "base_lattice_per_unit_depth" + # Force unit scales dx^3/dt^2, thus uniform fine force is twice root. + root_cd = 1 / (0.5 * r.SPEED**2 * 40) + fine_cd = 2 / (0.5 * r.SPEED**2 * 80) + assert root_cd == fine_cd + + +def test_pressure_local_frame_and_units(): + _, pair, _, _ = r.host_geometry(r.DEFAULT_SOURCE, h) + center = [ + (80.5 - pair.fine_origin[0]) / 0.5 - 0.5, + (81 - pair.fine_origin[1]) / 0.5 - 0.5, + ] + p = {"center": center, "radius": 40.0, "mean_speed": r.SPEED} + yy, xx = np.indices(pair.fine_shape) + rho = 1 + 1e-7 * xx + probes = h.pressure_probes(rho, np.ones(pair.fine_shape, bool), p) + assert np.isclose(probes["dp"], -80e-7 / (3 * r.SPEED**2)) + + +def test_initializer_coordinates_and_rest_ghosts(): + g = r.initial_state((166, 881), (0, 0), 1, root=True) + rho, ux, uy = r.macros(g) + assert np.max(abs(rho - 1)) < 1e-14 + assert np.max(abs(uy)) == 0 + assert np.max(abs(ux[[0, -1]])) == 0 + assert np.isclose(ux[81, 50], 6 * r.SPEED * (80.5 / 164) * (1 - 80.5 / 164)) + + +def test_axis_specific_sampling_phase(): + _, pair, _, _ = r.host_geometry(r.DEFAULT_SOURCE, h) + assert pair.fine_origin == (38.0, 22.0) + j = np.arange(pair.fine_shape[1]) + k = np.arange(pair.fine_shape[0]) + fine_x = (pair.fine_origin[0] + (j + 0.5) / 2 - 0.5) * 0.0025 + fine_y = (pair.fine_origin[1] + (k + 0.5) / 2 - 1) * 0.0025 + np.testing.assert_allclose(fine_x / 0.00125, j + 75.5, atol=1e-12) + np.testing.assert_allclose(fine_y / 0.00125, (k + 43) - 0.5, atol=1e-12) + np.testing.assert_allclose(fine_x - (j + 75) * 0.00125, 0.5 * 0.00125, atol=1e-15) + np.testing.assert_allclose(fine_y - ((k + 43) - 0.5) * 0.00125, 0, atol=1e-15) + + +def test_full_requires_completed_pilot(): + with unittest.TestCase().assertRaisesRegex(ValueError, "pilot-results"): + r.verify_pilot(None, r.verify_source(r.DEFAULT_SOURCE)) + + +def pilot_fixture(path): + """Synthetic contract fixture only; never a physical result or run receipt.""" + manifest = r.verify_source(r.DEFAULT_SOURCE) + (path / "plan.json").write_text( + json.dumps( + dict( + source_manifest=manifest, + bounds_root_edge=list(r.BOUNDS), + runner_sha256=r.runner_hash(), + ) + ) + ) + for case, diameter in (("uniform", 80), ("refined", 40)): + dest = path / case + dest.mkdir() + interval = round(0.02 * diameter / r.SPEED) + rows = [ + dict( + step=(i + 1) * interval, + tstar=(i + 1) * 0.02, + tstar_force=(i + 0.5) * 0.02, + cd=3.2, + cl=0.1, + dp=2.0, + mass=0.8, + rho_min=0.99, + rho_max=1.01, + max_mach=0.03, + flux_in=0.4, + flux_out=0.4, + rho_front=1.0, + rho_back=1.0, + ) + for i in range(250) + ] + (dest / "summary.json").write_text( + json.dumps( + dict( + case=case, + mode="pilot", + status="pilot_complete", + healthy=True, + steps=250 * interval, + history=rows, + startup_seconds=1.0, + wall_seconds=11.0, + step_seconds=9.0, + ) + ) + ) + (dest / "ready.json").write_text( + json.dumps( + dict( + case=case, + mode="pilot", + step=interval, + runner_sha256=r.runner_hash(), + all_runtime_and_observation_paths_exercised=True, + ) + ) + ) + return manifest + + +def test_pilot_validation_and_calibration(): + from tempfile import TemporaryDirectory + + with TemporaryDirectory() as tmp: + path = Path(tmp) + manifest = pilot_fixture(path) + assert r.verify_pilot(path, manifest) == 3606.0 + + +def test_reject_malformed_pilot_receipts(): + from tempfile import TemporaryDirectory + from copy import deepcopy + + with TemporaryDirectory() as tmp: + path = Path(tmp) + manifest = pilot_fixture(path) + receipt = path / "uniform/summary.json" + valid = json.loads(receipt.read_text()) + mutations = [ + lambda s: s.update(case="refined"), + lambda s: s.update(mode="full"), + lambda s: s.update(healthy="true"), + lambda s: s.update(steps=1), + lambda s: s["history"].pop(), + lambda s: s["history"][1].update(step=160), + lambda s: s["history"][-1].update(tstar=float("nan")), + lambda s: s["history"][0].update(cd=float("inf")), + lambda s: s["history"][0].update(max_mach=0.3), + lambda s: s["history"][0].update(tstar_force=0.02), + lambda s: s.update(wall_seconds=0), + ] + for mutate in mutations: + bad = deepcopy(valid) + mutate(bad) + receipt.write_text(json.dumps(bad)) + with unittest.TestCase().assertRaises(ValueError): + r.verify_pilot(path, manifest) + # Exact reviewer regression: previously accepted NaN endpoint. + receipt.write_text( + json.dumps( + dict( + status="pilot_complete", + healthy=True, + history=[{"tstar": float("nan")}], + ) + ) + ) + with unittest.TestCase().assertRaises(ValueError): + r.verify_pilot(path, manifest) + receipt.write_text(json.dumps(valid)) + (path / "uniform/ready.json").write_text("{}") + with unittest.TestCase().assertRaises(ValueError): + r.verify_pilot(path, manifest) + + +def test_strict_finite_json_and_deadline(): + from tempfile import TemporaryDirectory + + with TemporaryDirectory() as tmp: + p = Path(tmp) / "bad.json" + for data in ('{"x": NaN}', '{"x": Infinity}', '{"x": 1e999}'): + p.write_text(data) + with unittest.TestCase().assertRaises(ValueError): + r.read_finite_json(p) + with unittest.TestCase().assertRaises(TimeoutError): + r.check_deadline(0) + r.check_deadline(float("inf")) + + +def test_pressure_probes_actual_body_mask(): + from CelerisLab.lbm.descriptors import SOLID + from CelerisLab.lbm.refinement import ACTIVE + + _, pair, _, host = r.host_geometry(r.DEFAULT_SOURCE, h) + fluid = (pair.fine_kind == ACTIVE) & ( + (host.flags.reshape(pair.fine_shape) & SOLID) == 0 + ) + p = dict(center=host.objects[0].geometry.center, radius=40.0, mean_speed=r.SPEED) + assert h.pressure_probes(np.ones(pair.fine_shape), fluid, p)["dp"] == 0.0 + + +if __name__ == "__main__": + suite = unittest.TestSuite( + unittest.FunctionTestCase(fn) + for name, fn in sorted(globals().items()) + if name.startswith("test_") + ) + result = unittest.TextTestRunner(verbosity=2).run(suite) + raise SystemExit(not result.wasSuccessful()) diff --git a/tests/validation/data/refined_dfg/source_manifest.json b/tests/validation/data/refined_dfg/source_manifest.json new file mode 100644 index 0000000..5a65fd2 --- /dev/null +++ b/tests/validation/data/refined_dfg/source_manifest.json @@ -0,0 +1,110 @@ +{ + "commit": "cd3a657f9d3833e812c3cca173525f165b48710a", + "production_origin_commit": "47b05bc8abfc957fe876b0c63bdb386e3c6cbaa5", + "sha256": { + "src/CelerisLab/__init__.py": "21cc469bd6e970a1f9ff5970c8eea4602065e41fc881854759ee565ee55c10c6", + "src/CelerisLab/body/__init__.py": "e4e0e8f0769a94428d5b02fac790443997556cfc5515a5f1ed1595d2c694d69d", + "src/CelerisLab/body/action_smoother.py": "b1a23ee6bbd4850b172cb8738c5c3546ee465a54189e6af9d8044f023d4687dc", + "src/CelerisLab/body/coupling/__init__.py": "c931f7eb318d48408c3e2febcc6917c88c8dcc1697cbf4329b1f340d94c35763", + "src/CelerisLab/body/coupling/soa_packer.py": "8ec26fd55ef74801bbd5ab6bd3d0a7dabc750cbd498aaac2ca4c61a0108ad3c9", + "src/CelerisLab/body/ddf_patch.py": "60fcfb7441e73bf290bcff41ca3a83da3e7ddad765d7247464e8f02e3043b4e2", + "src/CelerisLab/body/geometry/__init__.py": "ada88980ab672959aca9aa3a4c4aade3219524b2ab68110cd2a2f41e3e0ca90a", + "src/CelerisLab/body/geometry/base.py": "cceca14361ca2f34d34a3e43730e1a79b8ca81be00f0dd62e884db32b479419f", + "src/CelerisLab/body/geometry/circle.py": "7605ca7c21103f1d801ac73ab8dd18dffbdac425b2dee7b8858b283628a9e436", + "src/CelerisLab/body/jets.py": "4104c24af4651f4bb4f8d8fd36fb6d88f41c6e4469822f36e22e74c9ede31dac", + "src/CelerisLab/body/manager.py": "01b653d25bd9b4333379c8bfd12383596b686d2029f2e91b1c91b618540ec6b3", + "src/CelerisLab/body/objects.py": "8e63cee8faa9cc6243bcff428be628dfad73a02a9a3b8ece7d78c0a63406a08d", + "src/CelerisLab/body/preprocess/__init__.py": "dfd9d58bc07f626ab23a70053e31b88a6d9165d2e4c76a4e9f48abe674d3d2ec", + "src/CelerisLab/body/preprocess/flag_overlay.py": "831fc9e9736b758847759a085ff133eb49ec91318ace4fda76f85960138ee300", + "src/CelerisLab/body/registry.py": "fbb7327780adbdc6da15ebdb799980daf32df22ee974af0893d7409f177bbb68", + "src/CelerisLab/body/sync_plan.py": "3a15bfe9ebeb91482dbfc68850761ef39c872cc12d2efc795babd0c01e392323", + "src/CelerisLab/common/__init__.py": "d4dbde78bc97fccc9eabac9abae8f0f15c112c10ff8e9b45bd9401d4723df1ef", + "src/CelerisLab/common/_types.py": "aad2ff1a9033d4fd68783ad7123afb47201ddc58785cb6dc0246947037072a5f", + "src/CelerisLab/common/checkpoint.py": "b2f960aecd5ac27a3a7130b809ffad4709c44d522edef0359ab952654625410b", + "src/CelerisLab/common/pathline.py": "a1b6efcc0a076d0f1dc8d49a595779426d7d5626b00072a59bac710ff02ad046", + "src/CelerisLab/common/preprocess.py": "7fc78a08857cf8e98519500b094a1d896af189b6b15658ff9b388ed51bed4515", + "src/CelerisLab/common/render.py": "0e055476257eb2bcd11e733766e52195de02f8fad228f570a3511a7c831bdc29", + "src/CelerisLab/common/streakline/__init__.py": "13ff1cfa416fdd721a08a659e76b5a5a99b8d3201add10a83b7b0bc654bc1eec", + "src/CelerisLab/common/streakline/_config.py": "a9185ef68e9c7b5d5a1dc0b04a28341e2ad91d13c145a7c106636a878cf0a605", + "src/CelerisLab/common/streakline/_integrate.py": "3fbd4368b2c8d0e7fda40dc65aac5df02b085008b4c68acd83cb7ae45d33a13b", + "src/CelerisLab/common/streakline/_render.py": "236cd637414ee318c4d66f9aed32110d8dc42a100d65fde06172f7cd0ba62c65", + "src/CelerisLab/common/streakline/_streakline.py": "1a3e08daa706b5f83cdeb337e9a945a36ddf7583c40af1d0bb9d825ac954df37", + "src/CelerisLab/config.py": "2f1bba7515f5bc08cb3a8d140bcafe7fa0689e01b3b06494697bb2956a71530e", + "src/CelerisLab/configs/CONFIG.md": "41cafdbe525366c08aa9d0416ea026ccf0026ba7123b8a9b65f968a80d928e3a", + "src/CelerisLab/configs/config_body.json": "055f261f0828ca1d4f2bb694bd15b72bbcabe31d552326208f28e37f761adab1", + "src/CelerisLab/configs/config_lbm.json": "17e4311f1f250c6ece59177a8b6ff856b69a016c7c10267e774b030a0a390865", + "src/CelerisLab/configs/config_lbm_flume_mrt_zouhe.json": "8f593321db05a563a42931603a228d61cbcba06dfbe113262d5cd0ad93c0fe79", + "src/CelerisLab/configs/config_lbm_flume_trt_regularized.json": "f4bba8e6a4339fc3c4eb14f6e2f57c4b22993c73703c773d853edfbe34736e56", + "src/CelerisLab/configs/config_lbm_three_cylinder_triangle.json": "eeb8e7149dbe4aa6b3a9562802735aa950750dfbe2d972a0727e788573afa817", + "src/CelerisLab/cuda/__init__.py": "0d01dff9f0399c8baba822aab5ca0122a135f011008ed8a499ebd9f7e5b06aea", + "src/CelerisLab/cuda/compiler_v2.py": "c4cc3ec4445c8d75e6a253dd87963877117a722bc73b9d4811806fc6045d08c1", + "src/CelerisLab/cuda/context.py": "013dd4ec6d01d2fbd11afdfb103eac6a273bd9f798d663d761b851eeb27885e9", + "src/CelerisLab/cuda/obs_layout.py": "5ffb2c3aeeaaaaeb64f2cd8febfbcf4ad9919e5414020ad753cfade0683b21c2", + "src/CelerisLab/lbm/__init__.py": "7f60fb3aaae81501b510d08c2fedff44fcb05657260eb10b6f5cc91e7204306b", + "src/CelerisLab/lbm/block_stepper.py": "00b313daddc6e63d8ab8423c9190ca7e92ad291de20a2f6b80eb53e3f2d85ec7", + "src/CelerisLab/lbm/blocks.py": "1e644c5719993f418185ad4a10b6950a30636fc912befd9eab41ff48b4cf85cd", + "src/CelerisLab/lbm/channel.py": "2d5b5ff43414cb35633ac4ce90ce2c3a3353950b25c168e1c2ab58d498dd8b0f", + "src/CelerisLab/lbm/curved_links.py": "b6957fa13d13f864dba4f7364e93a39905d062dbb3ee9b3f5191add81fff0732", + "src/CelerisLab/lbm/descriptors.py": "c9f025dfc14da43b989047b84cc22bd05d5cc92c9ec695c1da68fdfb788e3f2b", + "src/CelerisLab/lbm/equilibrium.py": "2fb2f2bb7ff326d23abce80fc2a4c0907b17740ea3e29fe6ae079c5044783826", + "src/CelerisLab/lbm/field.py": "0e727eacc3706add4d8d15c418e9b631050871f831698d3ac689dd5f2c701b9b", + "src/CelerisLab/lbm/hierarchy.py": "f234975074f6b6d289d573ff6ea6832feae0e27cb04f80b37fe81792b13a12e7", + "src/CelerisLab/lbm/initializers.py": "dacf9586963e862d6d0a98ef1a066036be6f91e9d423720661d8ef2c1fb2df6a", + "src/CelerisLab/lbm/kernels/blocks.cu": "e3b9606f861bfea8c4b0f34e8356bb230081f899c584f67be89178d889099b17", + "src/CelerisLab/lbm/kernels/boundary/bounce_back.cuh": "aba3748c564ed91b161d268ef8b737b7c7be96a823063d201f3da684dd845233", + "src/CelerisLab/lbm/kernels/boundary/curved_boundary.cuh": "aeb16bfd483514fb03e15a2fed437d673c3b8b2fe68ab83edf4244ef593878ea", + "src/CelerisLab/lbm/kernels/boundary/esopull_curved_helpers.cuh": "b9b474bc4edaff05681c5443235472d04072d0d775fc0eb9c5967faef52b0258", + "src/CelerisLab/lbm/kernels/boundary/ibm_kernels.cuh": "e6909b547dc8d64cdd8a76a90660243a8f5fee535caa161b52b4dd8b7c29934b", + "src/CelerisLab/lbm/kernels/boundary/inlet/channel_stabilized.cuh": "efe077afcbb6b731256a987c2fd79a1dbc470c902f87070d60421ea6cb13d4e2", + "src/CelerisLab/lbm/kernels/boundary/inlet/common.cuh": "59b38a89439878e00de7de16b9c7b03feff443fcc5a421593f7d2536dcbb7ab8", + "src/CelerisLab/lbm/kernels/boundary/inlet/equilibrium.cuh": "f0964c811b2c41a5aaa1c6d99c490b06e850c634eb8749f1f05ef61f93c2de1d", + "src/CelerisLab/lbm/kernels/boundary/inlet/regularized.cuh": "943d4a571818bacb89c70bccfaac202691121e3e385803ce6309123a0a03ddaf", + "src/CelerisLab/lbm/kernels/boundary/inlet/zou_he_local.cuh": "e73e66b853b1e4e0e30d8f4185d832aa414ce670a1813ce2b5459f2887f71f89", + "src/CelerisLab/lbm/kernels/boundary/inlet_outlet.cuh": "de8b213a6da42c2e5e9f8966199fc66ea4541b85144a90f3f5dbad3c6ab63b97", + "src/CelerisLab/lbm/kernels/boundary/outlet/pressure_neq.cuh": "c8b2adc66d12fe602076ddc9693c5c208711fcaddb58aa9a31401a269f51c4c4", + "src/CelerisLab/lbm/kernels/boundary/planar_velocity.cuh": "501b2759d0e2e4280cb6fdd127d66a50fa89f6de48fc2cae9a7eafb11ae90bed", + "src/CelerisLab/lbm/kernels/config/config_grid.h": "1db3e030ebdc116cd7f02362d1fd077daac55e8f0b4329f6285948c224d3ca5d", + "src/CelerisLab/lbm/kernels/config/config_method.h": "4374e04a130f90f97d17cc20c40d2f9ea0ed0dd8603817d07f215d02266e247c", + "src/CelerisLab/lbm/kernels/config/config_objects.h": "db793af052c756561a190a3ff2cb93ee5fe28462ee4f97c6af80882b4dc8b072", + "src/CelerisLab/lbm/kernels/config/config_obs.h": "c11e2249359145b5c4cbaaf0487506fa77c1f4b82c7c169ddc574aa0babcd727", + "src/CelerisLab/lbm/kernels/config/config_physics.h": "8ac0990b168e1058588676cc05c80df8634cfd671dbc92a89246c3c372d83a4e", + "src/CelerisLab/lbm/kernels/config.h": "7d1124125e0b354558297934460e64c760fb68cdd7bf6d01a1238bb2eaab8c37", + "src/CelerisLab/lbm/kernels/core/descriptors.cuh": "8440bc14857da9007982650f828a23608eec9f7cc4bf772e8c8dc1239cb45356", + "src/CelerisLab/lbm/kernels/core/flags.cuh": "0025705f24710330d9d9195f713879bfb650a85a3449baef0e196b831af104ab", + "src/CelerisLab/lbm/kernels/core/layout.cuh": "1fbd820021c89ac3efcde24145996cb0147769e607e498e61f308ed3e0b13706", + "src/CelerisLab/lbm/kernels/core/obs.cuh": "b742b77c9f61e1f020fe9922c3838ae9d855e17b381154958b58f10a22938922", + "src/CelerisLab/lbm/kernels/core/params.cuh": "6e18bf730d936d940e3d85452f0714a2eaae31bf41002f9fd520ec5b372878dd", + "src/CelerisLab/lbm/kernels/core/precision.cuh": "849df14b4cce50af52399e4b8f3e997196f3ae5f1163f68e28df3e6b2507a8c4", + "src/CelerisLab/lbm/kernels/kernel_v2.cu": "c7bd3fb705b77261ba2441c52a24a7e6a7865c165cba088a1b75201c4d29c408", + "src/CelerisLab/lbm/kernels/operators/collision_mrt.cuh": "aee444eb60c47a3fd518dca975d1af05425a0818a0ec6c3fb3f6150f60c2cb44", + "src/CelerisLab/lbm/kernels/operators/collision_srt.cuh": "6b7ff842f392c8f853d2d1a97cc69bfe50f21283f73e659404da07cf6c95f5e1", + "src/CelerisLab/lbm/kernels/operators/collision_trt.cuh": "9725b9840b9ad8b4ec7be61ff5d06d3151e86d1a87585c13f2fbf934d0d2c5f6", + "src/CelerisLab/lbm/kernels/operators/forcing_guo.cuh": "969a430c1af03f32e3bd0df95d805e64f56334717c44dd4ba208ba503062a24b", + "src/CelerisLab/lbm/kernels/operators/helpers.cuh": "0b3a4d95d7d19f5fec09e5423b8f4cc99f8a226b541967425967330da5a051ec", + "src/CelerisLab/lbm/kernels/operators/macro.cuh": "2be20eb2b7d7f5522daaad7f4b57b1ff8cdf883efe4a67846a96783b2fe846a3", + "src/CelerisLab/lbm/kernels/operators/turbulence_smag.cuh": "ff7f000297ee0817e255ad9fcb6a13b647e3ddd5b206518d35f38d2895a9f501", + "src/CelerisLab/lbm/kernels/refinement.cu": "129f3d91da93b4324f88a0da73a4da1e5390f72c14cce574ee63306712da9895", + "src/CelerisLab/lbm/kernels/step/aux_kernels.cu": "578a5b20a157241ee678547059abb1372dd5a589ad06bb2c15b987eb01ca1013", + "src/CelerisLab/lbm/kernels/step/block_channel.cu": "befbf6d757e532a59c166d41f920d16d71557cd3a574478e7aa83da737c48cbd", + "src/CelerisLab/lbm/kernels/step/curved_double.cu": "bfa7a1df1e939020314d49755045f7a85fe4a19cc1c68dbe62653c66c5e7c2b3", + "src/CelerisLab/lbm/kernels/step/esopull_macro.cu": "c73116020ca217ede1e4de27c8269c31f9d91e33433ac8252a21e9b14187c50f", + "src/CelerisLab/lbm/kernels/step/esopull_sync.cu": "c0ff3e72900c3ab89a05dd03684e96d3b740cf1672d565727fc2af165933d087", + "src/CelerisLab/lbm/kernels/step/init_flow.cu": "305be3bb12fc0d039acec1701baabcf0e8a35e683c49c34f2af8ca9867d2b29c", + "src/CelerisLab/lbm/kernels/step/one_step_double.cu": "f9f8242b38e8e0cf8425b75e2b0215bc15363c62eaa0fb43c761b8e0931bfdc9", + "src/CelerisLab/lbm/kernels/step/one_step_esopull.cu": "aa223b2744ae6c9125673df311a96ce588032c1ba58bebb11fd33ecb990f8ed4", + "src/CelerisLab/lbm/kernels/streaming/esopull_semantic_helpers.cuh": "dd5a9091628b079efd97e53fd49e52d09dcee88529235db77a5648e02937af6e", + "src/CelerisLab/lbm/kernels/streaming/esopull_single_buffer.cuh": "77ed4ab9a0ebc5ff62a91a7c39066ef395391ae36e9040a69fb2bd8b1626dc79", + "src/CelerisLab/lbm/kernels/streaming/pull_double_buffer.cuh": "b65d741ed78b7b62f99531ba0c0ce1936663d8dc9640e7d7d16c707f1baf6549", + "src/CelerisLab/lbm/population_storage.py": "a6a86c29c52c892edfe4802fb773c387212e3f4cc86708d80e154986ddf88044", + "src/CelerisLab/lbm/raw_exchange.py": "8694202466eba848563e298656353b0f5db94e21a01fa05f51aad6dc2b9fd311", + "src/CelerisLab/lbm/refinement.py": "08c2a0557276a98cbcaa459d9bb7a9735d8bbbf64a6114c57af2fda61481a687", + "src/CelerisLab/lbm/refinement_stepper.py": "f2bb34fa90805ad21c3e7cf9f9a67d88d3ae3e249b4517d3898567b4bf9e364d", + "src/CelerisLab/lbm/refinement_tables.py": "fdd344024b89df3a83174dc810ca5e32db5e19eb8b1609f8f73b0821f56546b6", + "src/CelerisLab/lbm/stepper.py": "0dcb7cd7ae87e33e3dc67a38b6aa43042fbe673e37d9ad165342ee1d695933ce", + "src/CelerisLab/lbm/synchronized_remap.py": "7400fc4beaffd4649d6a9c524d38ce5fd6fa7eeb32eb5332e4030a4b6648f197", + "src/CelerisLab/simulation.py": "0a03904dc6de88667d7af699005a289adf8b09f0976cb0d099b041f4830a8f4c", + "tests/validation/re100_helpers.py": "289c6e01075ff361573f862b90cfd61e5980c9a8efcbc78aef69979b4004fd60", + "tests/validation/periodic.py": "a2212c0e2705126327a5019e9cafe64eeca398da1b75061ff23630c9843697a6", + "tests/validation/data/dfg_re100/reference_cycle.json": "073d1bda1629697efebad95667bfe80598e217d708dbffd7d5f6c6999379c55f" + } +} \ No newline at end of file diff --git a/tests/validation/run_refined_dfg.py b/tests/validation/run_refined_dfg.py new file mode 100644 index 0000000..cc4acff --- /dev/null +++ b/tests/validation/run_refined_dfg.py @@ -0,0 +1,603 @@ +"""Frozen matched SRT DFG comparison. GPU execution requires explicit root admission.""" + +import argparse +import csv +import hashlib +import json +import math +import os +from pathlib import Path +import sys +import time +from tempfile import TemporaryDirectory + +import numpy as np + +HERE = Path(__file__).resolve().parent +DEFAULT_SOURCE = HERE.parents[1] +BOUNDS = (40, 24, 320, 140) +BASE_BLOCK = (2, 881) # Exact 166 x 881 domain; no geometric padding. +FINE_BLOCK = (8, 8) # 240 x 568 including four fine halo cells. +SPEED = 0.01 + + +def helpers(source): + sys.path.insert(0, str(source / "src")) + sys.path.insert(0, str(source / "tests/validation")) + import CelerisLab + + if ( + Path(CelerisLab.__file__).resolve() + != source.resolve() / "src/CelerisLab/__init__.py" + ): + raise RuntimeError("Wrong CelerisLab import") + import re100_helpers as h + import periodic + + if h.ROOT.resolve() != source.resolve(): + raise RuntimeError("Helper module loaded from wrong source") + return h, periodic + + +def config(h, diameter): + cfg, body, p = h.make_config(diameter, SPEED, "SRT", "neq_extrap", "zou_he_local") + cfg["method"]["outlet"]["srt_neq_damp"] = 0.5 + return cfg, body, p + + +def plan(h): + _, _, c = config(h, 40) + _, _, f = config(h, 80) + x0, y0, x1, y1 = BOUNDS + fine_active = 4 * (x1 - x0) * (y1 - y0) + coarse_uncovered = c["nx"] * c["ny"] - fine_active // 4 + return dict( + status="prepared_not_run", + root_gpu_gates_required=True, + root=c, + uniform=f, + bounds_root_edge=list(BOUNDS), + physical_patch_m=[ + (x0 - 0.5) * c["dx"], + (y0 - 1) * c["dx"], + (x1 - 0.5) * c["dx"], + (y1 - 1) * c["dx"], + ], + body_root_edge=[80.5, 81.0], + radius_root=20.0, + root_block=list(BASE_BLOCK), + fine_block=list(FINE_BLOCK), + root_wall_edges=[1.0, 165.0], + root_inlet_outlet_edges=[0.5, 880.5], + physical_domain_m=[2.2, 0.41], + physical_center_m=[0.2, 0.2], + physical_diameter_m=0.1, + physical_mean_speed_m_s=1.0, + physical_viscosity_m2_s=0.001, + finest_dx_m=f["dx"], + finest_dt_s=f["dt"], + approximate_leaf_cell_ratio=(coarse_uncovered + fine_active) + / (f["nx"] * f["ny"]), + approximate_collision_work_ratio=(coarse_uncovered + 2 * fine_active) + / (2 * f["nx"] * f["ny"]), + cost_caveat="Geometry-only estimates include ghosts/solids; exclude block halos, transfers, observation and launch overhead. Not measured speedup.", + sampling_offset="For origin (38,22), fine physical x=(j+75.5)*dx80 versus uniform i*dx80: half-fine-cell x offset. Fine physical y=(k+42.5)*dx80 equals uniform (l-.5)*dx80 at l=k+43: y samples align. Body, wall and source planes coincide.", + outlet="SRT pressure NEQ, all nine fneq damped by .5, no backflow clamp; noncanonical DFG outlet", + modes={"pilot_tstar": 5.0, "full_tstar": 300.0}, + observation_dense_tstar=0.02, + ) + + +def macros(g): + rho = np.sum(g, axis=0, dtype=np.float64) + ux = (g[1] - g[2] + g[5] - g[6] + g[7] - g[8]) / rho + uy = (g[3] - g[4] + g[5] - g[6] - g[7] + g[8]) / rho + return rho, ux, uy + + +def leaf_mass(arrays, dx, masks): + """Finite public leaves only; covered coarse, solid and halos are excluded.""" + if any(not np.isfinite(g[:, mask]).all() for g, mask in zip(arrays, masks)): + raise ValueError("Nonfinite physical leaf population") + return sum( + float(np.sum(np.sum(g, axis=0, dtype=np.float64)[mask], dtype=np.float64)) + * (dx / 2**level) ** 2 + for level, (g, mask) in enumerate(zip(arrays, masks)) + ) + + +def initial_state(shape, origin, spacing, height=164.0, speed=SPEED, root=False): + """Same physical parabolic initial condition as uniform InitTubeFlow_v2.""" + y = origin[1] + spacing * (np.arange(shape[0]) + 0.5) + u = 6 * speed * ((y - 1) / height) * (1 - (y - 1) / height) + velocity = np.broadcast_to(u[:, None], shape).copy() + if root: + velocity[[0, -1], :] = 0 + velocity[:, [0, -1]] = 0 + cx = np.array([0, 1, -1, 0, 0, 1, -1, 1, -1])[:, None, None] + w = np.array([4 / 9] + [1 / 9] * 4 + [1 / 36] * 4)[:, None, None] + cu = cx * velocity + return w * (1 + 3 * cu + 4.5 * cu * cu - 1.5 * velocity * velocity) + + +def host_geometry(source, h): + """Production donor/transfer and complete cut-link checks, without CUDA.""" + from CelerisLab.config import load_lbm_config + from CelerisLab.lbm.refinement import FixedPatch2D, prepare_fine_boundary + from CelerisLab.lbm.blocks import BlockLayout2D + from CelerisLab.lbm.channel import ( + validate_channel_config, + validate_channel_transfers, + ) + from CelerisLab.lbm.refinement_tables import blocked_base_tables + from CelerisLab.body.geometry.circle import CircleGeometry + from CelerisLab.body.objects import SimObject + + cfgdict, _, _ = config(h, 40) + with TemporaryDirectory() as tmp: + p = Path(tmp) / "cfg.json" + p.write_text(json.dumps(cfgdict)) + cfg = load_lbm_config(str(p)) + pair = FixedPatch2D.from_config(cfg, BOUNDS) + layout = BlockLayout2D.build(cfg.ny, cfg.nx, *BASE_BLOCK) + validate_channel_config(cfg, [BASE_BLOCK, FINE_BLOCK]) + targets = blocked_base_tables(pair, layout, np.argsort(pair.coalesce_coarse))[1] + validate_channel_transfers(pair, layout, targets) + BlockLayout2D.build(*pair.fine_shape, *FINE_BLOCK, periodic=(False, False)) + body = SimObject( + 0, geometry=CircleGeometry(80.5, 81.0, 20.0), center=(80.5, 81.0), radius=20.0 + ) + host = prepare_fine_boundary(cfg, pair, (body,), _shifted_blocked=True) + return cfg, pair, body, host + + +def verify_source(source): + expected = json.loads((HERE / "data/refined_dfg/source_manifest.json").read_text()) + for name, digest in expected["sha256"].items(): + if hashlib.sha256((source / name).read_bytes()).hexdigest() != digest: + raise RuntimeError("Source mismatch: " + name) + return expected + + +def run_case(a, h, periodic, case): + from CelerisLab.simulation import Simulation + from CelerisLab.cuda.context import CudaContext + from CelerisLab.lbm.refinement_stepper import FixedPatchStepper + from CelerisLab.lbm.descriptors import FLUID + from contextlib import ExitStack + + refined = case == "refined" + cfgdict, bodydict, p = config(h, 40 if refined else 80) + cap = 5.0 if a.mode == "pilot" else 300.0 + max_steps = math.ceil(cap * p["diameter"] / SPEED) + rows = [] + healthy = True + step_seconds = 0.0 + started = time.monotonic() + check_deadline(a.deadline) + dest = a.output / case + dest.mkdir() + with ExitStack() as stack: + if refined: + cfg, pair, body, host = host_geometry(a.source, h) + context = stack.enter_context(CudaContext(0)) + backend = stack.enter_context( + FixedPatchStepper( + cfg, + BOUNDS, + context=context, + coarse_boundary="channel", + base_block=BASE_BLOCK, + fine_block=FINE_BLOCK, + fine_objects=(body,), + observe_raw_exchange=True, + ) + ) + from CelerisLab.lbm.channel import channel_masks + from CelerisLab.lbm.refinement import ACTIVE + from CelerisLab.lbm.descriptors import SOLID + + masks = [ + channel_masks(pair.coarse_shape)[0] & ~pair.covered, + (pair.fine_kind == ACTIVE) + & ((host.flags.reshape(pair.fine_shape) & SOLID) == 0), + ] + backend.initialize( + initial_state(pair.coarse_shape, (0.0, 0.0), 1.0, root=True), + initial_state(pair.fine_shape, pair.fine_origin, 0.5), + ) + probe_p = dict( + p, + center=[ + (80.5 - pair.fine_origin[0]) / 0.5 - 0.5, + (81.0 - pair.fine_origin[1]) / 0.5 - 0.5, + ], + radius=40.0, + ) + else: + tmp = Path(stack.enter_context(TemporaryDirectory())) + (tmp / "lbm.json").write_text(json.dumps(cfgdict)) + (tmp / "body.json").write_text(json.dumps(bodydict)) + sim = stack.enter_context( + Simulation(str(tmp / "lbm.json"), str(tmp / "body.json"), device_id=0) + ) + sim.initialize() + probe_p = p + fluid = (sim.field.flag.reshape(p["ny"], p["nx"]) & FLUID) != 0 + startup = time.monotonic() - started + check_deadline(a.deadline) + begin = 0 + with (dest / "history.csv").open("w") as out: + writer = None + while begin < max_steps: + check_deadline(a.deadline) + width = ( + 0.02 + if a.mode == "pilot" or begin * SPEED / p["diameter"] >= 240 + else 0.5 + ) + n = min(max(1, round(width * p["diameter"] / SPEED)), max_steps - begin) + impulse = np.zeros(2) + tick = time.monotonic() + if refined: + for _ in range(n): + check_deadline(a.deadline) + backend.step(1) + exchange = backend.read_raw_exchange( + 0, level=1, window="last_base_step" + ) + if ( + exchange["units"] != "base_lattice_per_unit_depth" + or exchange["sample_count"] != 2 + or exchange["time_end"] - exchange["time_start"] != 1 + ): + raise RuntimeError("Unexpected force observation contract") + impulse += exchange["total"]["impulse"] + check_deadline(a.deadline) + else: + done = 0 + while done < n: + check_deadline(a.deadline) + batch = min(100, n - done) + sim.run(batch) + impulse += sim.read_force(0).astype(np.float64) * batch + done += batch + check_deadline(a.deadline) + step_seconds += time.monotonic() - tick + if refined: + arrays = backend.download() + fields = [macros(g) for g in arrays] + root_rho, root_ux, _ = fields[0] + rho, ux, uy = fields[1] + fluid = masks[1] + if any( + not np.isfinite(g[:, mask]).all() + for g, mask in zip(arrays, masks) + ): + healthy = False + break + mass = leaf_mass(arrays, p["dx"], masks) + vals = [ + tuple(v[mask] for v in field) + for mask, field in zip(masks, fields) + ] + else: + m = sim.get_macroscopic() + rho, ux, uy = [ + np.asarray(m[k], np.float64).reshape(p["ny"], p["nx"]) + for k in ("rho", "ux", "uy") + ] + root_rho, root_ux = rho, ux + mass = float(rho[fluid].sum()) * p["dx"] ** 2 + vals = [(rho[fluid], ux[fluid], uy[fluid])] + healthy = ( + all(np.isfinite(v).all() for vs in vals for v in vs) + and np.isfinite(impulse).all() + ) + if not healthy: + break + den = SPEED**2 * p["diameter"] / 2 + row = dict( + step=begin + n, + tstar=(begin + n) * SPEED / p["diameter"], + tstar_force=(begin + n / 2) * SPEED / p["diameter"], + cd=float(impulse[0] / n / den), + cl=float(impulse[1] / n / den), + mass=mass, + rho_min=float(min(vs[0].min() for vs in vals)), + rho_max=float(max(vs[0].max() for vs in vals)), + max_mach=float( + max(np.sqrt(3 * (vs[1] ** 2 + vs[2] ** 2)).max() for vs in vals) + ), + flux_in=float(np.sum(root_rho[1:-1, 2] * root_ux[1:-1, 2])) + * p["dx"] + / SPEED, + flux_out=float(np.sum(root_rho[1:-1, -3] * root_ux[1:-1, -3])) + * p["dx"] + / SPEED, + **h.pressure_probes(rho, fluid, probe_p), + ) + check_deadline(a.deadline) + ensure_finite(row) + rows.append(row) + if writer is None: + writer = csv.DictWriter(out, fieldnames=row) + writer.writeheader() + writer.writerow(row) + out.flush() + healthy = ( + row["max_mach"] <= 0.2 + and row["rho_min"] >= 0.5 + and row["rho_max"] <= 1.5 + ) + begin += n + if not healthy: + break + if len(rows) == 1: + ready = dict( + case=case, + mode=a.mode, + step=begin, + tstar=row["tstar"], + runner_sha256=runner_hash(), + all_runtime_and_observation_paths_exercised=True, + ) + (dest / "ready.json").write_text( + json.dumps(ready, indent=2, allow_nan=False) + ) + print("READY " + json.dumps(ready, allow_nan=False), flush=True) + result = dict( + case=case, + mode=a.mode, + history=rows, + healthy=bool(healthy), + steps=begin, + startup_seconds=startup, + step_seconds=step_seconds, + wall_seconds=time.monotonic() - started, + force="raw stationary-body link momentum reaction, per-unit-depth; average over entire observation interval; normalized in each base lattice", + mass_units="rho_reference * m^2 per unit depth", + flux_units="rho_reference * m^2/s per unit depth", + pressure="postcollision end-sample; physical pressure difference / (rho_reference * 1 m/s squared)", + suite_validated=False, + **periodic.assess(rows), + ) + if not healthy: + result["periodic"] = False + reference = json.loads( + (a.source / "tests/validation/data/dfg_re100/reference_cycle.json").read_text() + ) + if result["cycles"]: + result["reference_comparison"] = periodic.compare( + result["cycles"][-1], reference + ) + result["status"] = ( + "invalid_state" + if not healthy + else ( + "pilot_complete" + if a.mode == "pilot" + else ("periodic" if result["periodic"] else "incomplete_periodicity") + ) + ) + (dest / "summary.json").write_text(json.dumps(result, indent=2, allow_nan=False)) + return result + + +def runner_hash(): + return hashlib.sha256(Path(__file__).read_bytes()).hexdigest() + + +def ensure_finite(value): + if isinstance(value, dict): + for item in value.values(): + ensure_finite(item) + elif isinstance(value, list): + for item in value: + ensure_finite(item) + elif isinstance(value, (float, int)) and not isinstance(value, bool): + if not math.isfinite(value): + raise ValueError("Nonfinite JSON value") + + +def read_finite_json(path): + def reject(value): + raise ValueError("Nonfinite JSON constant: " + value) + + value = json.loads(path.read_text(), parse_constant=reject) + ensure_finite(value) + return value + + +def check_deadline(deadline): + if time.monotonic() >= deadline: + raise TimeoutError( + "Campaign wall-time ceiling reached; no continuation or acceptance" + ) + + +def verify_pilot(path, manifest): + if path is None: + raise ValueError( + "Full comparison requires --pilot-results from a completed matched pilot" + ) + previous = read_finite_json(path / "plan.json") + if previous["source_manifest"] != manifest or previous["bounds_root_edge"] != list( + BOUNDS + ): + raise ValueError("Pilot source or frozen geometry differs") + if previous.get("runner_sha256") != runner_hash(): + raise ValueError("Pilot runner differs") + calibrated_seconds = 0.0 + for case, diameter in (("uniform", 80), ("refined", 40)): + summary = read_finite_json(path / case / "summary.json") + if ( + summary.get("case") != case + or summary.get("mode") != "pilot" + or summary.get("status") != "pilot_complete" + or summary.get("healthy") is not True + ): + raise ValueError("Mislabeled or unhealthy pilot receipt") + rows = summary.get("history") + interval = round(0.02 * diameter / SPEED) + total = round(5 * diameter / SPEED) + if ( + not isinstance(rows, list) + or len(rows) != 250 + or summary.get("steps") != total + ): + raise ValueError("Incomplete pilot history or step count") + for i, row in enumerate(rows): + step = (i + 1) * interval + if type(row.get("step")) is not int or row["step"] != step: + raise ValueError("Nonmonotone or incomplete pilot step history") + for key, expected in ( + ("tstar", step * SPEED / diameter), + ("tstar_force", (step - interval / 2) * SPEED / diameter), + ): + value = row.get(key) + if type(value) not in (float, int) or not math.isclose( + value, expected, rel_tol=0, abs_tol=1e-12 + ): + raise ValueError("Incorrect pilot observation time") + for key in ( + "cd", + "cl", + "dp", + "mass", + "rho_min", + "rho_max", + "max_mach", + "flux_in", + "flux_out", + "rho_front", + "rho_back", + ): + if type(row.get(key)) not in (int, float): + raise ValueError("Missing numeric pilot observation") + if not ( + 0 < row["mass"] + and 0.5 <= row["rho_min"] <= row["rho_max"] <= 1.5 + and 0 <= row["max_mach"] <= 0.2 + ): + raise ValueError("Unhealthy pilot observations") + ready = read_finite_json(path / case / "ready.json") + if ( + ready.get("case") != case + or ready.get("mode") != "pilot" + or ready.get("step") != interval + or ready.get("runner_sha256") != runner_hash() + or ready.get("all_runtime_and_observation_paths_exercised") is not True + ): + raise ValueError("Invalid pilot READY receipt") + startup = summary.get("startup_seconds") + wall = summary.get("wall_seconds") + stepping = summary.get("step_seconds") + if ( + any(type(v) not in (float, int) for v in (startup, wall, stepping)) + or not 0 <= startup < wall + or not 0 < stepping <= wall + ): + raise ValueError("Invalid pilot calibration timings") + # Full physical duration is 60 times pilot; dense pilot sampling makes + # this conservative for the sparser early full observations. Not a guarantee. + calibrated_seconds += startup + 60 * (wall - startup) + return 3 * calibrated_seconds + + +def main(): + parser = argparse.ArgumentParser(description=__doc__) + parser.add_argument("--source", type=Path, default=DEFAULT_SOURCE) + parser.add_argument("--output", type=Path, required=True) + parser.add_argument("--mode", choices=["plan", "pilot", "full"], default="plan") + parser.add_argument("--gpu-gates-approved", action="store_true") + parser.add_argument("--pilot-results", type=Path) + parser.add_argument( + "--max-wall-seconds", + type=float, + help="Required explicit campaign wall-time ceiling for GPU execution", + ) + a = parser.parse_args() + a.source = a.source.resolve() + manifest = verify_source(a.source) + h, periodic = helpers(a.source) + calibration = verify_pilot(a.pilot_results, manifest) if a.mode == "full" else None + if a.mode != "plan" and ( + a.max_wall_seconds is None + or not math.isfinite(a.max_wall_seconds) + or a.max_wall_seconds <= 0 + ): + raise ValueError("GPU execution requires finite positive --max-wall-seconds") + budget = ( + min(a.max_wall_seconds, calibration) + if calibration is not None + else a.max_wall_seconds + ) + if a.output.exists() and any(a.output.iterdir()): + raise ValueError("Output must be empty/new") + a.output.mkdir(parents=True, exist_ok=True) + specification = plan(h) + specification.update( + source_manifest=manifest, + gates=periodic.GATES, + reference_gates=periodic.REFERENCE_GATES, + runner_sha256=runner_hash(), + wall_time_budget_seconds=budget, + full_pilot_calibration_seconds=calibration, + ) + host_geometry(a.source, h) + (a.output / "plan.json").write_text(json.dumps(specification, indent=2)) + if a.mode == "plan": + return 0 + if not a.gpu_gates_approved or not os.environ.get("CUDA_VISIBLE_DEVICES"): + raise RuntimeError( + "GPU execution blocked: root GPU gates approval and explicit device selection required" + ) + results = [] + a.deadline = time.monotonic() + budget + for case in ("uniform", "refined"): + try: + result = run_case(a, h, periodic, case) + except TimeoutError as exc: + failure = dict( + status="wall_time_limit", + case=case, + matched_case_pass=False, + error=str(exc), + ) + (a.output / "comparison.json").write_text( + json.dumps(failure, indent=2, allow_nan=False) + ) + print(json.dumps(failure, allow_nan=False), flush=True) + return 3 + results.append(result) + if not result["healthy"]: + break + comparison = dict( + status="pilot_only" if a.mode == "pilot" else "incomplete", + matched_case_pass=False, + ) + if len(results) == 2 and all(r["healthy"] and r["periodic"] for r in results): + comparison.update( + status="periodic_comparison", + refined_vs_uniform=periodic.compare( + results[1]["cycles"][-1], results[0]["cycles"][-1] + ), + ) + comparison["matched_case_pass"] = comparison["refined_vs_uniform"][ + "reference_agreement" + ] + (a.output / "comparison.json").write_text( + json.dumps(comparison, indent=2, allow_nan=False) + ) + return ( + 0 + if len(results) == 2 + and all(r["healthy"] for r in results) + and (a.mode == "pilot" or comparison["matched_case_pass"]) + else 2 + ) + + +if __name__ == "__main__": + raise SystemExit(main())