diff --git a/experiments/REDUCED_SURFACE_STUDY_2026-09-18.md b/experiments/REDUCED_SURFACE_STUDY_2026-09-18.md index 3a6daee..ddd0feb 100644 --- a/experiments/REDUCED_SURFACE_STUDY_2026-09-18.md +++ b/experiments/REDUCED_SURFACE_STUDY_2026-09-18.md @@ -1,5 +1,7 @@ # One-hour reduced smooth surface study — 2026-09-18 +Outcome: independent smooth surface parameters remove the immediate full-step obstruction and reduce cold linear solve times by 2–5 times in this screen. They do not remove the original shape-force imbalance. A three-step seed chain reaches P(h) virial errors below 7e-7 while retaining substantial excluded force and imperfect interior fields. + Authorized wall-clock window: 19:30:15–20:30:15 UTC. Branch: research/reduced-smooth-surface-20260918, parent 6ceda61. @@ -49,6 +51,13 @@ Experiment-only headers: excluded-force decomposition, raw trials and restoration invariants. - `audit_reduced_surface.py`: independent artifact algebra and corruption tests. +For physical state factors S and residual factors R, the shape maps obey + + S B = Qx Tx, R^-1 B = Qf Tf, + Qf^T (R Fphysical) = Tf^-T B^T Fphysical. + +The separate state/test QR maps therefore preserve physical dual pairing despite nonuniform normalization and distinct shape row/value offsets. The area whitening measures scalar radial amplitude, not the norm of the generated three-dimensional displacement. The latter is measured separately over the stellar volume. Changing the area whitening changes coordinates without changing the raw harmonic span. + The diagnostic executable adds experiment-only flags --reduced-surface-probe and --surface-degree (0 means full-space control); no production-facing API is modified. Degrees are runtime experimental choices while operator/inverse @@ -59,6 +68,108 @@ necessarily checked at runtime. Existing physical observers are reused. Evidence root: `../../diagnostics/reduced_surface_2026-09-18` relative to this checkout. The preceding study remains on its two committed research branches. -## Results +## Results: late-state comparison -Pending the bounded serial campaign and its independent artifact audit. +At the common six-update n=1.5, J=0 Q3 endpoint, the normalized original residual is 1.4971971e-5. The unrestricted control reproduces the previous fresh direction. Its unit trial is rejected at pressure preparation (thermodynamic domain); alpha=1/4 and 1/16 prepare but increase the original norm. Only alpha=1/1024 passes among the sampled values. Intermediate steps were not tested. + +| Shape space | Independent parameters | Krylov iterations | Cold solve seconds | Original relative linear defect | Reduced relative linear defect | Original residual after unit trial | +|---|---:|---:|---:|---:|---:|---:| +| Full Q3 trace | 866 | 61 | 32.184 | 0.02767 | 0.02767 | Pressure preparation rejected | +| L2 harmonics | 9 | 12 | 6.447 | 0.15121 | 0.01719 | 2.263903e-6 | +| L4 harmonics | 25 | 21 | 11.249 | 0.16187 | 0.01608 | 2.432960e-6 | +| L8 harmonics | 81 | 29 | 15.745 | 0.19207 | 0.01409 | 2.964119e-6 | + +All reduced directions pass their independent 3% linear gate, and all four trial alphas prepare and decrease both norms. None of the reduced directions passes the original full Newton gate. Original unit-step residual decreases are 84.88%, 83.75%, and 80.20%. Thus the benefit is measurable in the original residual, but it comes from solving a restricted set of shape equations. These directions do not meet the requested 3% accuracy for the original Newton system; they could be interpreted as much looser inexact full-system steps. + +These are one-shot cold linear solve timings, not end-to-end Newton timings or repeated performance statistics. The remaining non-shape unknowns and equations stay active: full dimension 178073, L4 dimension 177232. The full operator and full inverse are still applied inside the reduced wrappers. The dimension change is less than 0.5%; fewer Krylov iterations provide most of the measured speedup. + +### Curvature and physical motion + +| Shape space | Volume displacement direction L2 norm (q14) | Unit-trial original Taylor remainder | +|---|---:|---:| +| Full | 3.96028e-3 | Inadmissible; alpha=1/4 remainder1.85956e-4 | +| L2 | 4.86622e-7 | 1.84872e-10 | +| L4 | 9.14977e-5 | 2.21090e-7 | +| L8 | 1.46119e-4 | 7.19540e-7 | + +L4 cuts physical displacement motion by 43.28 times relative to the full direction. L2 is nearly a field-only relaxation; admitting higher harmonics introduces substantially more compensated shape/field motion and nonlinear curvature. The degree dependence is not monotone improvement of a physical solution. This supports the approximation-space compatibility hypothesis without proving that every omitted mode is unphysical or identifying a continuum gauge. + +After the L4 unit trial, excluded shape forces account for 98.38% of the remaining original squared residual. A single trial does not establish a new residual floor. Because q0 is retained, this late-state test also retains its existing roughly 1.1% nodal radius spread; restricting subsequent increments cannot remove the general fine-scale part of that inherited geometry. + +### Physical measurements + +The table below compares the common base with unit trials, all using independent q18 evaluation. Pressure is evaluated through enthalpy, P(h). These energy- and force-based virial diagnostics must not be confused with the density-pressure diagnostic, which remains undefined where the interpolated density is negative. + +| Measurement | Base | L2 unit | L4 unit | L8 unit | +|---|---:|---:|---:|---:| +| P(h) energy virial error | 4.22590e-6 | 8.66263e-7 | 8.11913e-7 | 8.15910e-7 | +| P(h) force virial error | 2.71367e-6 | 5.94930e-7 | 6.41398e-7 | 6.24187e-7 | +| Mass error relative to 1 | 3.41685e-6 | 8.81237e-10 | 2.39038e-8 | 5.67517e-8 | +| EOS density mismatch RMS | 1.82986e-4 | 1.82925e-4 | 1.82990e-4 | 1.83010e-4 | +| Hydrostatic gradient RMS | 1.06811e-2 | 1.06631e-2 | 1.07065e-2 | 1.06864e-2 | +| Minimum sampled density | -7.73185e-4 | -7.24527e-4 | -7.27628e-4 | -7.33640e-4 | +| Negative-density mass fraction | 2.89683e-5 | 2.89242e-5 | 2.89458e-5 | 2.89452e-5 | + +Both P(h) virial balances are comfortably below 1e-5. Their q14-to-q18 differences are about 1.83e-9. Local EOS and hydrostatic errors remain essentially unchanged; the negative density has not been repaired. Negative-density mass estimates are quadrature-sensitive (about 2.23e-5 at q14 versus2.89e-5 at q18), so their magnitude is diagnostic rather than a precision estimate. All quoted gradient errors are the existing interior observers in code units, including the broken DG-potential gradient comparison where applicable. + +The full-control physical.csv sample is alpha=1/4, its first prepared trial, which fails the norm screen. It is not a selected full-space update and must not be compared as an equivalent accepted unit step. + +## Seed continuation + +The planned analytic-seed L4 probe was followed by two externally replayed physical unit candidates, an adaptive extension recorded before launching the chain. Each predecessor passed preparation, the reduced independent linear gate, and both conservative norm-decrease screens. All three observed unit trials satisfy the same conditions. This is a cold-solve experimental reduced-Newton chain, not the production accepted/warm trajectory. The first physical shape is retained, and each subsequent base is its preceding saved candidate. + +| Unit steps from seed | Original norm | Retained norm | Excluded shape norm | P(h) energy virial q18 | Hydrostatic gradient RMS q18 | Nodal radius spread / mean | +|---|---:|---:|---:|---:|---:|---:| +| 0 | 2.736592e-4 | 2.590732e-4 | 8.815010e-5 | 2.410694e-5 | 1.528098e-3 | 0.003424% | +| 1 | 8.624827e-5 | 9.387845e-6 | 8.573583e-5 | 2.180785e-6 | 5.724618e-3 | 0.081116% | +| 2 | 8.473526e-5 | 3.147803e-6 | 8.467677e-5 | 8.660498e-7 | 7.177698e-3 | 0.152605% | +| 3 | 8.465191e-5 | 7.357283e-8 | 8.465188e-5 | 6.769653e-7 | 7.255981e-3 | 0.154699% | + +The retained norm falls by about 3521 times. The original norm falls by only 3.23 times overall; after the first step, improvements are 1.75% and 0.0984%. At the final trial, excluded shape force accounts for 99.999924% of the original squared norm. These three updates show a pronounced separation between progress in the retained equations and the remaining force; they do not establish an asymptotic residual floor or convergence of the retained problem to arbitrary tolerance. + +Linear solves take 18/20/18 iterations and 10.574/11.341/9.594 seconds. Reduced relative defects are 2.940%, 2.122%, 2.323%; original relative defects are 31.45%, 98.18%, 99.90%. All full-system 3% gates fail. + +The final P(h) force virial is 5.769064e-7, mass error 3.79142e-9, sampled minimum density -4.07472e-4, and negative-density mass fraction 2.89249e-5. EOS density mismatch RMS scaled improves from 1.72399e-4 at seed to 1.25532e-4, but the hydrostatic-gradient error grows 4.75 times. Radius spread includes the reference geometry's initial interpolation error; this nodal statistic is not a full surface norm. The restricted chain produces much less shape distortion than the older unrestricted six-update state, but does not eliminate nonrotating drift or establish pointwise correctness. + +This is the decisive limitation of a permanently truncated surface space: good virial balance and a small retained norm coexist with substantial original shape force and imperfect interior fields. The late-state trial benefits from fine-scale balance already developed by the unrestricted trajectory. Smooth parameters relieve the immediate curvature/step obstruction; they do not by themselves establish a compatible approximation of all original equations. + + +![Original and retained residual observations](../../../diagnostics/reduced_surface_2026-09-18/reduced_surface_results.png) + +## Interpretation and next design + +The evidence supports separating independent surface parameters from the interpolation used to carry volume displacement. An internal typed discretization prescription can select a nested smooth surface space while preserving the physics-facing API. State injection, physical dual restriction, normalization, and excluded-force monitoring should be distinct interfaces. A reduced inverse can initially wrap the existing full inverse; later optimizations should preserve the same contracts. + +A fixed small harmonic space is not yet justified for production. It changes which shape forces are enforced, and smaller reduced residuals do not certify the full equations. The next controlled study should test whether the excluded mechanical imbalance decreases under mesh or field enrichment at a spherical state, while measuring interior EOS/hydrostatic error. That can distinguish approximation error from an inconsistent force/operator coupling before a low-degree space is promoted to production. Degree and mesh sensitivity, the pressure/EOS representation, and compatibility of the volume extension with global radial dilation remain separate questions. This study makes no rotation or n>1.5 convergence claim and changes no production field order or volume lift. + +## Verification and reproducibility + +Implementation commit 21cdbda22486e408adc52b65b25413f9380c7c30, on research/reduced-smooth-surface-20260918, parent 6ceda61. The diagnostic includes 19 C++ algebra checks (maximum error 2.22e-16) and 22 independent Python algebra/corruption tests. Actual-case audits independently reconstruct the raw nodal harmonics, quadrature-defined scalar trace Gram matrix, both QR maps, physical dual identity, full/reduced defects and norm decomposition, all trial remainder/norm classifications, and exact saved-state restoration. Individual artifact audits do not independently reimplement the FEM physics observers. + +The observation flags use norm(Ftrial)<=(1-1e-4 alpha)norm(Fbase). This is a deliberately conservative norm-decrease screen, stricter than the production baseline squared-merit criterion, not an exact replay of production AcceptStep. Physical measurements cover the base and the first prepared trial only. No probe commits a production accepted state, and no SOC or positivity repair is enabled. + +All builds, FEM cases and large artifact audits were serial. Raw commands, mesh/input/output and executable hashes, compiler settings, runtime and CPU details are in each campaign provenance and run JSON. The source checkout stays on the isolated branch; the original checkout and its pre-existing edits are untouched. Full regression tests were not rerun because production sources are unchanged; verification is scoped to the new diagnostic contracts, actual cases and snapshot chain. + + +All seven case audits pass: 2,098 saved-artifact checks in total. The snapshot-chain audit passes 107 checks with all three individual reviews bound. Both physical-state and residual handoffs are bitwise identical; normalized-state round trips differ by at most 1.39e-17 per entry. Layout, state/residual factors and raw harmonic basis hashes match exactly. Shape-base increments remain in the original raw span to relative errors below 9.4e-16. Twelve small chain-auditor tests include corruption and JSON-serialization regressions. + +The chain auditor initially failed while serializing a NumPy boolean; it also emitted warnings during its small BLAS matrix product. The checker was corrected to emit native JSON scalar types and use an explicitly checked elementwise contraction with finite-value guards. The corrected audit passes without warnings. No FEM state, recorded observation, producer binary or numerical source changed in response to this checker issue. + +### Archived campaigns + +Evidence directory (relative to the source checkout): +../../diagnostics/reduced_surface_2026-09-18/ + +| Campaign | Cases | Case audit | +|---|---|---| +| core_v1 | late_full_control, late_L4 | core_audit_v1.json | +| orders_v1 | late_L2, late_L8 | order_audit_v1.json | +| seed_v1 | seed_L4 | seed_audit_v1.json | +| seed_chain_v1 | seed_L4_step2 | seed_step2_audit_v1.json | +| seed_chain_v2 | seed_L4_step3 | seed_step3_audit_v1.json | + +seed_chain_audit_v1.json binds the three seed stages. TEST_PLAN.md records the original scope and the conditional extension. late_summary.json and seed_summary.json collect the displayed values. plot_results.py regenerates the figure from the immutable case CSVs. COMPLETED_STUDY.json records final source and evidence hashes. + +To reproduce, build plateau_diagnostic from producer commit 21cdbda using the compiler/MFEM configuration in core_v1/provenance.json. Replay each recorded run command, substituting a fresh output directory. All cases use the same saved sandbox.smesh and one MPI rank. Late cases bind the prior Q3 accepted_state.txt by hash. Replay the seed stages in order using physical_candidate_a0 from the previous stage; the case manifests and input-state hashes specify these links. The runner requires a clean tracked source checkout and records the executable and all case output hashes. + +No order changes were made during this study. Optional L6, rotational models, higher indices, mesh/enrichment sweeps and production integration were outside this bounded experiment. All seven FEM probes finished by 20:09:58 UTC, before the 20:20 launch cutoff and 20:30:15 study deadline. diff --git a/experiments/audit_reduced_surface_chain.py b/experiments/audit_reduced_surface_chain.py new file mode 100644 index 0000000..024d753 --- /dev/null +++ b/experiments/audit_reduced_surface_chain.py @@ -0,0 +1,377 @@ +#!/usr/bin/env python3 +"""Check up to three ordered frozen reduced-surface probes as an affine chain. + +This checks saved handoffs and the common physical shape span. Individual FEM +artifact audits remain separate; --case-audit binds their completed JSON reviews +without rerunning their larger matrix checks. Passing this audit establishes +neither an unrestricted Newton solve nor a physically valid stellar endpoint. +""" +import argparse +import csv +import hashlib +import json +import math +import os +from pathlib import Path +import sys + +for _name in ("OPENBLAS_NUM_THREADS", "OMP_NUM_THREADS", "MKL_NUM_THREADS", "VECLIB_MAXIMUM_THREADS"): + os.environ[_name] = "1" +import numpy as np + +PREFIX = "reduced_surface_" +EPS = np.finfo(float).eps + + +def digest(path): + result = hashlib.sha256() + with Path(path).open("rb") as stream: + for chunk in iter(lambda: stream.read(1024 * 1024), b""): + result.update(chunk) + return result.hexdigest() + + +def rows(path): + with Path(path).open(newline="") as stream: + reader = csv.DictReader(stream) + if not reader.fieldnames or len(set(reader.fieldnames)) != len(reader.fieldnames): + raise ValueError(f"Missing or duplicate CSV columns: {path}") + result = list(reader) + if any(None in row or any(value is None for value in row.values()) for row in result): + raise ValueError(f"Malformed CSV data: {path}") + return result + + +def finite(value): + result = float(value) + if not math.isfinite(result): + raise ValueError(f"Nonfinite value: {value}") + return result + + +def integer(value): + result = finite(value) + if result != int(result): + raise ValueError(f"Nonintegral value: {value}") + return int(result) + + +def load(path, matrix=False): + with Path(path).open() as stream: + shape = tuple(int(value) for value in stream.readline().split()) + values = np.loadtxt(stream, ndmin=1) + if len(shape) != (2 if matrix else 1) or min(shape) <= 0 or values.ndim != 1 or values.size != math.prod(shape): + raise ValueError(f"Invalid saved array dimensions: {path}") + if not np.isfinite(values).all(): + raise ValueError(f"Nonfinite saved array: {path}") + return values.reshape(shape) + + +def norm(values): + # Stable scaling keeps an overflowing squared norm from becoming an + # apparently permissive infinite threshold in the handoff checks. + scale = float(np.max(np.abs(values), initial=0.)) + result = 0. if scale == 0. else scale * math.sqrt(float(np.sum((values / scale) ** 2))) + return finite(result) + + +def compare(actual, expected, atol, rtol): + if actual.shape != expected.shape: + return {"passed": False, "reason": "shape_mismatch"} + difference = actual - expected + bound = atol + rtol * np.maximum(np.abs(actual), np.abs(expected)) + maximum = float(np.max(np.abs(difference), initial=0.)) + scale = max(norm(actual), norm(expected)) + return {"passed": bool(np.all(np.abs(difference) <= bound)), + "bitwise_equal": actual.tobytes() == expected.tobytes(), + "maximum_absolute_difference": maximum, "difference_norm": norm(difference), + "relative_difference_norm": norm(difference) / scale if scale else 0., + "absolute_tolerance": atol, "relative_tolerance": rtol} + + +def span_check(basis, increment, reference, atol, rtol): + coefficients, _, rank, singular = np.linalg.lstsq(basis, increment, rcond=None) + if not np.isfinite(coefficients).all() or not np.isfinite(singular).all(): + raise ValueError("Nonfinite least-squares coefficients or singular values") + # This small validation contraction needs no BLAS dispatch. Explicit + # products also make a real overflow distinguishable from platform BLAS + # floating-status warnings observed despite finite saved input vectors. + with np.errstate(over="raise", invalid="raise", divide="raise"): + reconstructed = np.sum(basis * coefficients[None, :], axis=1) + if not np.isfinite(reconstructed).all(): + raise ValueError("Nonfinite reconstructed shape increment") + residual = increment - reconstructed + error, size = norm(residual), norm(increment) + # Subtracting two full shape bases has a small absolute rounding floor; + # report it separately rather than inflating the relative span tolerance. + subtraction_floor = float(64 * EPS * max(norm(reference), norm(reference + increment), 1e-300)) + threshold = float(atol * math.sqrt(increment.size) + subtraction_floor + rtol * size) + return {"passed": bool(int(rank) == basis.shape[1] and error <= threshold), + "rank": int(rank), "requested_rank": basis.shape[1], "increment_norm": size, + "projection_error_norm": error, "relative_projection_error": error / size if size else 0., + "projection_error_threshold": threshold, "subtraction_roundoff_allowance": subtraction_floor, + "largest_singular_value": float(singular[0]), "smallest_singular_value": float(singular[-1]), + "coefficients": coefficients.tolist()} + + +def inspect(directories, atol=1e-14, rtol=5e-10, span_rtol=5e-10, reviews=()): + failures, checks, cases, links = [], 0, [], [] + + def check(label, condition): + nonlocal checks + checks += 1 + if not condition: + failures.append(label) + + individual = {} + review_hashes = {} + try: + folders = [Path(path).resolve() for path in directories] + if not 1 <= len(folders) <= 3 or len(set(folders)) != len(folders): + raise ValueError("Provide one to three distinct case directories in chain order") + for review_path in reviews: + review_path = Path(review_path).resolve() + review = json.loads(review_path.read_text()) + review_hashes[str(review_path)] = digest(review_path) + for case in review.get("cases", [review]): + if "directory" not in case: + raise ValueError("Individual review lacks its case directory") + key = str(Path(case["directory"]).resolve()) + if key in individual: + raise ValueError(f"Duplicate individual case review: {key}") + individual[key] = case + + for index, folder in enumerate(folders): + used = {} + + def file(name, extension=".txt"): + path = folder / (PREFIX + name + extension) + used[path.name] = digest(path) + return path + + summary_rows = rows(file("summary", ".csv")) + summary = {row["metric"]: finite(row["value"]) for row in summary_rows} + if len(summary) != len(summary_rows): + raise ValueError("Duplicate summary metrics") + layout = rows(file("layout", ".csv")) + shape_rows = [row for row in layout if row["kind"] == "value" and row["block"] == "surface_deformation"] + if len(shape_rows) != 1: + raise ValueError("Missing or duplicate physical shape value block") + offset, size = (integer(shape_rows[0][key]) for key in ("offset", "size")) + data = {name: load(file(name)) for name in + ("physical_base", "base", "residual", "state_factors", "residual_factors", "reduced_residual", "reduced_defect")} + basis = load(file("basis_raw"), matrix=True) + n = data["base"].size + check(f"case{index}:shape_layout", 0 <= offset <= n - size and size > 0 and basis.shape[0] == size) + check(f"case{index}:array_dimensions", all(data[key].shape == (n,) for key in + ("physical_base", "residual", "state_factors", "residual_factors"))) + check(f"case{index}:positive_factors", bool(np.all(data["state_factors"] > 0) and np.all(data["residual_factors"] > 0))) + roundtrip = compare(data["base"], data["state_factors"] * data["physical_base"], atol, rtol) + check(f"case{index}:physical_normalized_roundtrip", roundtrip["passed"]) + check(f"case{index}:reduced_model", summary["identity_control"] == 0 and + summary["shape_value_offset"] == offset and summary["shape_size"] == size and + summary["mode_count"] == basis.shape[1]) + reduced_norm, defect_norm = norm(data["reduced_residual"]), norm(data["reduced_defect"]) + relative, absolute = summary["linear_relative_tolerance"], summary["linear_absolute_tolerance"] + check(f"case{index}:linear_tolerance", 0 < relative < 1 and absolute >= 0 and reduced_norm > 0) + threshold = max(relative * reduced_norm, absolute) + reduced_gate = defect_norm <= threshold + check(f"case{index}:reported_reduced_gate", summary["reduced_verified_gate"] in (0, 1) and + bool(summary["reduced_verified_gate"]) == reduced_gate) + + trials = rows(file("trials", ".csv")) + unit_rows = [row for row in trials if integer(row["sample"]) == 0] + if len(unit_rows) != 1 or finite(unit_rows[0]["alpha"]) != 1.: + raise ValueError("Expected exactly one alpha=1 sample zero") + unit = unit_rows[0] + unit_prepared = integer(unit["prepared"]) + check(f"case{index}:unit_prepared_boolean", unit_prepared in (0, 1)) + candidate = None + merit = {"prepared": bool(unit_prepared)} + if unit_prepared: + candidate = {name: load(file(name + "_a0")) for name in + ("physical_candidate", "candidate", "residual", "reduced_residual")} + full_pass = norm(candidate["residual"]) <= (1. - 1e-4) * norm(data["residual"]) + reduced_pass = norm(candidate["reduced_residual"]) <= (1. - 1e-4) * reduced_norm + for label, passed in (("full", full_pass), ("reduced", reduced_pass)): + flag = integer(unit[label + "_merit_pass"]) + check(f"case{index}:reported_{label}_unit_merit", flag in (0, 1) and bool(flag) == passed) + merit.update(full_merit_pass=full_pass, reduced_merit_pass=reduced_pass, + full_norm=norm(candidate["residual"]), reduced_norm=norm(candidate["reduced_residual"])) + else: + check(f"case{index}:unprepared_merit_flags", integer(unit["full_merit_pass"]) == integer(unit["reduced_merit_pass"]) == 0) + qualifies = bool(reduced_gate and unit_prepared and merit.get("full_merit_pass") and merit.get("reduced_merit_pass")) + # Only predecessor probes must offer a successful handoff. A final + # failed or rejected unit trial honestly terminates a valid prefix. + if index < len(folders) - 1: + check(f"case{index}:qualifies_for_next_base", qualifies) + + run_path, provenance_path = folder.parent / (folder.name + "_run.json"), folder.parent / "provenance.json" + run = json.loads(run_path.read_text()) + provenance = json.loads(provenance_path.read_text()) + check(f"case{index}:completed_producer", bool(run.get("finished_utc")) and run.get("exit_code") == 0 and not run.get("timed_out", False)) + for name, observed in used.items(): + check(f"case{index}:producer_hash:{name}", run.get("output_hashes", {}).get(name) == observed) + binding = individual.get(str(folder)) + bound = False + if reviews: + check(f"case{index}:individual_review_present", binding is not None) + if binding is not None: + source = binding.get("observations", {}).get("provenance", {}) + bound = bool(binding.get("artifact_pass") and source.get("bound") and + source.get("run_sha256") == digest(run_path) and + source.get("provenance_sha256") == digest(provenance_path)) + check(f"case{index}:individual_review_bound", bound) + record = {"directory": str(folder), "input_hashes": used, "run_sha256": digest(run_path), + "provenance_sha256": digest(provenance_path), "individual_audit_bound": bound, + "normalization_roundtrip": roundtrip, "reduced_verified_gate": reduced_gate, + "reduced_defect_norm": defect_norm, "reduced_gate_threshold": threshold, + "unit_trial": merit, "qualifies_for_continuation": qualifies, + "physical_shape_span": None} + if index == 0: + first_basis, first_shape = basis, data["physical_base"][offset:offset + size].copy() + shared_hashes = {name: used[PREFIX + name + extension] for name, extension in + (("layout", ".csv"), ("state_factors", ".txt"), ("residual_factors", ".txt"), ("basis_raw", ".txt"))} + first_provenance = provenance + else: + for name, expected in shared_hashes.items(): + extension = ".csv" if name == "layout" else ".txt" + check(f"case{index}:identical_{name}", used[PREFIX + name + extension] == expected) + for key in ("binary_sha256", "mesh_sha256"): + check(f"case{index}:identical_{key}", bool(provenance.get(key)) and provenance.get(key) == first_provenance.get(key)) + if previous_candidate is None: + raise ValueError("A chain predecessor has no prepared unit candidate") + comparisons = {label: compare(data[current], previous_candidate[prior], atol, rtol) for label, current, prior in + (("physical", "physical_base", "physical_candidate"), ("normalized", "base", "candidate"), + ("residual", "residual", "residual"))} + for label, observation in comparisons.items(): + check(f"handoff{index}:{label}", observation["passed"]) + links.append({"from": str(folders[index - 1]), "to": str(folder), **comparisons}) + observation = span_check(first_basis, data["physical_base"][offset:offset + size] - first_shape, + first_shape, atol, span_rtol) + check(f"case{index}:same_affine_shape_span", observation["passed"]) + record["physical_shape_span"] = observation + cases.append(record) + previous_candidate = candidate + except (OSError, ValueError, KeyError, TypeError, IndexError, FloatingPointError, np.linalg.LinAlgError) as error: + check("complete_consistent_chain_inputs:" + str(error), False) + passed = not failures + all_bound = bool(cases) and len(cases) == len(directories) and all(case["individual_audit_bound"] for case in cases) + return {"schema_version": 1, "passed": passed, "checks_count": checks, "failures": failures, + "individual_audits_all_bound": all_bound, "complete_artifact_validation": passed and all_bound, + "case_audit_hashes": review_hashes, "cases": cases, "handoffs": links, + "tolerances": {"handoff_absolute": atol, "handoff_relative": rtol, "shape_span_relative": span_rtol}, + "scope": "Frozen saved-vector affine chain only; individual artifact audits are separate. No new FEM evaluation, unrestricted convergence, admissibility certificate, or physical validation is claimed."} + + +def self_test(): + # A good endpoint difference can be small in total norm yet outside the + # allowed shape span. The orthogonal perturbation must still be caught. + basis = np.array([[1., 0.], [0., 1.], [1., 1.]]) + reference = np.array([.1, .2, .3]) + increment = basis @ np.array([.03, -.01]) + span = span_check(basis, increment, reference, 1e-14, 5e-10) + assert span["passed"] + assert json.loads(json.dumps(span, allow_nan=False))["passed"] is True + assert not span_check(basis, increment + 1e-5 * np.array([1., 1., -1.]), reference, 1e-14, 5e-10)["passed"] + assert not span_check(np.column_stack((basis[:, 0], basis[:, 0])), increment, reference, 1e-14, 5e-10)["passed"] + assert compare(reference, reference.copy(), 1e-14, 5e-10)["bitwise_equal"] + roundtrip = compare(reference + 2e-15, reference, 1e-14, 5e-10) + assert roundtrip["passed"] and not roundtrip["bitwise_equal"] + assert not compare(reference + 1e-5, reference, 1e-14, 5e-10)["passed"] + # Two tiny producer folders exercise the handoff logic itself, including + # the distinction between an honest stopped prefix and a falsely extended + # chain. These are not substitutes for independent per-case FEM audits. + import tempfile + with tempfile.TemporaryDirectory() as temporary: + root = Path(temporary) + (root / "provenance.json").write_text(json.dumps({"binary_sha256": "synthetic", "mesh_sha256": "synthetic"})) + + def csv_save(folder, name, columns, values): + with (folder / (PREFIX + name + ".csv")).open("w", newline="") as stream: + writer = csv.writer(stream) + writer.writerow(columns) + writer.writerows(values) + + def save(folder, name, array): + array = np.asarray(array) + with (folder / (PREFIX + name + ".txt")).open("w") as stream: + stream.write(" ".join(str(value) for value in array.shape) + "\n") + np.savetxt(stream, array.ravel(), fmt="%.17g") + + def seal(folder): + record = {"finished_utc": "synthetic_complete", "exit_code": 0, + "output_hashes": {path.name: digest(path) for path in folder.iterdir()}} + (root / (folder.name + "_run.json")).write_text(json.dumps(record)) + + first, second = root / "first", root / "second" + for index, folder in enumerate((first, second)): + folder.mkdir() + base = np.array([1., .1 * index, .1 * index, 1.]) + candidate = base + np.array([0., .1, .1, 0.]) + residual = np.full(4, .5 ** index) + reduced = np.full(3, .5 ** index) + for name, array in {"physical_base": base, "base": base, "residual": residual, + "state_factors": np.ones(4), "residual_factors": np.ones(4), + "reduced_residual": reduced, "reduced_defect": .01 * reduced, + "basis_raw": np.ones((2, 1)), "physical_candidate_a0": candidate, + "candidate_a0": candidate, "residual_a0": .5 * residual, + "reduced_residual_a0": .5 * reduced}.items(): + save(folder, name, array) + csv_save(folder, "summary", ("metric", "value"), + {"identity_control": 0, "shape_value_offset": 1, "shape_size": 2, "mode_count": 1, + "linear_relative_tolerance": .03, "linear_absolute_tolerance": 0., "reduced_verified_gate": 1}.items()) + csv_save(folder, "layout", ("kind", "block", "offset", "size"), + (("value", "head", 0, 1), ("value", "surface_deformation", 1, 2), ("value", "tail", 3, 1))) + csv_save(folder, "trials", ("sample", "alpha", "prepared", "full_merit_pass", "reduced_merit_pass"), + ((0, 1, 1, 1, 1),)) + seal(folder) + result = inspect((first, second)) + assert result["passed"], result["failures"] + assert json.loads(json.dumps(result, allow_nan=False))["passed"] is True + assert not result["complete_artifact_validation"] + original = (first / (PREFIX + "reduced_defect.txt")).read_bytes() + save(first, "reduced_defect", np.ones(3)) + seal(first) + result = inspect((first, second)) + assert not result["passed"] and "case0:qualifies_for_next_base" in result["failures"] + (first / (PREFIX + "reduced_defect.txt")).write_bytes(original) + seal(first) + save(second, "physical_base", np.array([1., .10001, .1, 1.])) + seal(second) + result = inspect((first, second)) + assert not result["passed"] and "handoff1:physical" in result["failures"] + print("12 tiny affine-span, JSON, handoff, complete-prefix and corruption checks passed") + + +def main(): + parser = argparse.ArgumentParser(description=__doc__) + parser.add_argument("directories", nargs="*", type=Path) + parser.add_argument("--output", type=Path) + parser.add_argument("--case-audit", type=Path, action="append", default=[]) + parser.add_argument("--absolute-tolerance", type=float, default=1e-14) + parser.add_argument("--relative-tolerance", type=float, default=5e-10) + parser.add_argument("--span-relative-tolerance", type=float, default=5e-10) + parser.add_argument("--self-test", action="store_true") + args = parser.parse_args() + if args.self_test: + self_test() + return 0 + if not 1 <= len(args.directories) <= 3 or args.output is None: + parser.error("Provide one to three ordered case directories and --output") + tolerances = (args.absolute_tolerance, args.relative_tolerance, args.span_relative_tolerance) + if any(not math.isfinite(value) or value < 0 for value in tolerances): + parser.error("Tolerances must be finite and nonnegative") + output = args.output.resolve() + if output.exists() or any(output.is_relative_to(path.resolve()) for path in args.directories): + parser.error("Use a fresh output outside immutable case directories") + result = inspect(args.directories, *tolerances, reviews=args.case_audit) + result["checker_sha256"] = digest(__file__) + output.parent.mkdir(parents=True, exist_ok=True) + output.write_text(json.dumps(result, indent=2, allow_nan=False) + "\n") + print(f"{'PASS' if result['passed'] else 'FAIL'}: {result['checks_count']} chain checks; {output}") + return 0 if result["passed"] else 1 + + +if __name__ == "__main__": + sys.exit(main())