From a7d646d04979692f780cdaf00d6c1ec7310323fc Mon Sep 17 00:00:00 2001 From: Emily Boudreaux Date: Thu, 17 Sep 2026 09:44:04 -0400 Subject: [PATCH] Check frozen numerical strategies without conflating incomplete evidence with success --- experiments/check_numerical_strategies.py | 315 ++++++++++++++++++ .../check_numerical_strategies_README.md | 72 ++++ .../test_check_numerical_strategies.py | 170 ++++++++++ 3 files changed, 557 insertions(+) create mode 100644 experiments/check_numerical_strategies.py create mode 100644 experiments/check_numerical_strategies_README.md create mode 100644 experiments/test_check_numerical_strategies.py diff --git a/experiments/check_numerical_strategies.py b/experiments/check_numerical_strategies.py new file mode 100644 index 0000000..6cf8a91 --- /dev/null +++ b/experiments/check_numerical_strategies.py @@ -0,0 +1,315 @@ +#!/usr/bin/env python3 +"""Check frozen numerical-strategy probes, without certifying a stellar model. + +Exit 0: requested evidence passes. Exit 3: an invariant or expectation fails. +Exit 2: evidence is incomplete (including a correctly reported failed solve). +Use check_plateau_physics.py separately for an end-to-end equilibrium endpoint. +""" + +import argparse +import csv +import json +import math +from pathlib import Path +import sys + + +EPSILON = sys.float_info.epsilon +SIGN_STATUSES = {"certified_nonnegative", "unresolved", "negative", "non_finite"} + + +def read_rows(path, required): + with path.open(newline="") as stream: + reader = csv.DictReader(stream) + if reader.fieldnames is None or len(set(reader.fieldnames)) != len(reader.fieldnames): + raise ValueError(f"Missing or duplicate columns in {path.name}") + if not set(required).issubset(reader.fieldnames): + raise ValueError(f"Missing required columns in {path.name}: {sorted(set(required) - set(reader.fieldnames))}") + result = [] + for row in reader: + if None in row or any(value is None for value in row.values()): + raise ValueError(f"Malformed CSV row in {path.name}") + result.append(row) + return result + + +def metrics(path): + result = {} + for row in read_rows(path, ("metric", "value")): + if not row["metric"] or row["metric"] in result: + raise ValueError(f"Empty or duplicate metric in {path.name}") + result[row["metric"]] = float(row["value"]) + return result + + +def number(row, key): + return float(row[key]) + + +def flag(row, key): + if row[key] not in ("0", "1"): + raise ValueError(f"Invalid boolean {key}: {row[key]!r}") + return row[key] == "1" + + +def nonnegative(value): + return math.isfinite(value) and value >= 0.0 + + +def roundoff(*values): + return 512.0 * EPSILON * (1.0 + sum(abs(value) for value in values if math.isfinite(value))) + + +def clean_json(value): + if isinstance(value, float) and not math.isfinite(value): + return None + if isinstance(value, dict): + return {key: clean_json(item) for key, item in value.items()} + if isinstance(value, list): + return [clean_json(item) for item in value] + return value + + +class Assessment: + def __init__(self, folder, restoration_tolerance): + self.folder = folder + self.restoration_tolerance = restoration_tolerance + self.checks = [] + self.outcomes = {} + + def check(self, name, condition, **details): + self.checks.append({"name": name, "status": "pass" if condition else "fail", **details}) + + def incomplete(self, name, **details): + self.checks.append({"name": name, "status": "incomplete", **details}) + + def restore(self, values, accepted, restored): + for key in accepted: + self.check(key, values[key] == 0.0, value=values[key], required=0.0) + for key in restored: + self.check(key, nonnegative(values[key]) and values[key] <= self.restoration_tolerance, + value=values[key], tolerance=self.restoration_tolerance) + + def polynomial(self, expect_negative): + certificates = read_rows(self.folder / "polynomial_certificates.csv", ("state", "status", "witness_element", "x", "y", "z")) + by_state = {row["state"]: row for row in certificates} + if len(by_state) != len(certificates): + raise ValueError("Duplicate polynomial certificate state") + for name in ("frozen_base", "context_seed_anchor", "original_full_candidate"): + if name not in by_state: + raise ValueError(f"Missing polynomial certificate: {name}") + if any(row["status"] not in SIGN_STATUSES for row in certificates): + raise ValueError("Unknown polynomial certificate status") + repairs = read_rows(self.folder / "polynomial_repair.csv", ("status", "minimum_safe_fraction", "base_status")) + if len(repairs) != 1: + raise ValueError("Polynomial repair summary must contain one row") + repair = repairs[0] + trials = read_rows(self.folder / "polynomial_repair_trials.csv", ("fraction", "alpha", "certificate", "prepared", "merit_slope")) + base_status = by_state["frozen_base"]["status"] + self.outcomes["polynomial"] = {"base": base_status, "anchor": by_state["context_seed_anchor"]["status"], + "repair": repair["status"], "trials": len(trials)} + self.check("repair_base_matches_certificate", repair["base_status"] == base_status) + if expect_negative: + self.check("expected_negative_base_detected", base_status == "negative", observed=base_status) + self.check("negative_base_repair_rejected", repair["status"] == "invalid_base", observed=repair["status"]) + if base_status != "certified_nonnegative": + self.check("uncertified_base_is_rejected", repair["status"] == "invalid_base") + self.check("uncertified_base_has_no_repaired_trials", all(number(row, "fraction") == 0.0 for row in trials)) + if base_status == "unresolved": + self.incomplete("polynomial_base_classification", reason="Subdivision budget did not resolve the base sign") + if base_status == "negative": + witness = by_state["frozen_base"] + self.check("negative_base_has_local_witness", int(witness["witness_element"]) >= 0 and + all(0.0 <= number(witness, axis) <= 1.0 for axis in ("x", "y", "z"))) + if repair["status"] == "feasible": + fraction = number(repair, "minimum_safe_fraction") + self.check("repair_fraction_in_unit_interval", math.isfinite(fraction) and 0.0 <= fraction <= 1.0, + fraction=fraction) + repaired = [row for row in trials if number(row, "fraction") > 0.0] + if fraction > 0.0 and not repaired: + self.incomplete("certified_repair_trials", reason="No repaired candidates were evaluated") + for index, row in enumerate(repaired): + self.check(f"repaired_polynomial_{index}_certified", row["certificate"] == "certified_nonnegative") + self.outcomes["polynomial"]["prepared_repaired_trials"] = sum(flag(row, "prepared") for row in repaired) + self.outcomes["polynomial"]["descent_repaired_trials"] = sum(number(row, "merit_slope") < 0.0 for row in repaired) + elif repair["status"] != "invalid_base": + self.incomplete("polynomial_repair_available", outcome=repair["status"]) + restoration = metrics(self.folder / "polynomial_restoration.csv") + self.restore(restoration, + ("accepted_physical_difference", "accepted_normalized_difference", "accepted_residual_difference", "restored_correction_difference"), + ("restored_residual_difference",)) + + def gravity(self): + summary = metrics(self.folder / "numerical_strategy_summary.csv") + self.restore(summary, ("accepted_state_change", "accepted_residual_change", "accepted_physical_change"), + ("maximum_original_restoration_error",)) + inner = read_rows(self.folder / "numerical_strategy_inner.csv", ("kind", "alpha", "index", "status", "converged", "right_hand_side_norm", "true_residual_norm", "relative_true_residual_norm", "threshold")) + lookup = {} + for index, row in enumerate(inner): + key = (row["kind"], number(row, "alpha"), int(row["index"])) + if key in lookup: + raise ValueError(f"Duplicate inner solve record: {key}") + lookup[key] = row + converged = flag(row, "converged") + self.check(f"inner_{index}_status_consistent", (int(row["status"]) == 0) == converged) + if converged: + rhs, residual, threshold = (number(row, key) for key in ("right_hand_side_norm", "true_residual_norm", "threshold")) + self.check(f"inner_{index}_true_defect_gate", all(nonnegative(x) for x in (rhs, residual, threshold)) and residual <= math.nextafter(threshold, math.inf), + true_residual=residual, threshold=threshold) + expected_relative = residual / rhs if rhs > 0 else (0.0 if residual == 0 else math.inf) + self.check(f"inner_{index}_relative_defect", math.isclose(number(row, "relative_true_residual_norm"), expected_relative, rel_tol=1e-12, abs_tol=1e-30)) + else: + self.incomplete(f"inner_{index}_solve_available", kind=row["kind"], alpha=number(row, "alpha"), + outcome="numerical_solve_failure", solve_status=int(row["status"])) + self.outcomes["gravity_inner"] = {"recorded_calls": len(inner), "failed_calls": sum(not flag(row, "converged") for row in inner)} + self.curvature(summary, lookup) + self.reduced(summary, lookup) + + def curvature(self, summary, lookup): + rows = read_rows(self.folder / "numerical_strategy_curvature.csv", ("alpha", "prepared", "inner_converged", "corrected_prepared", "raw_gravity_remainder", "corrected_gravity_remainder", "gravity_correction_norm", "correction_over_alpha_squared", "retained_state_change")) + if len(rows) != 4 or sorted(number(row, "alpha") for row in rows) != [.001, .01, .1, 1.]: + self.incomplete("curvature_probe_coverage", expected_alphas=[1, .1, .01, .001], recorded_rows=len(rows)) + reliable = [] + evaluated = 0 + for row in rows: + alpha = number(row, "alpha") + if not math.isfinite(alpha) or alpha <= 0: + raise ValueError("A curvature step must be finite and positive") + prepared, solved, corrected = (flag(row, key) for key in ("prepared", "inner_converged", "corrected_prepared")) + self.check(f"curvature_{alpha:g}_stages_consistent", (not solved or prepared) and (not corrected or solved)) + if not (prepared and solved and corrected): + self.incomplete(f"curvature_{alpha:g}_corrected_trial_available", raw_prepared=prepared, + inner_converged=solved, corrected_prepared=corrected) + continue + evaluated += 1 + solve = lookup[("curvature", alpha, 0)] + self.check(f"curvature_{alpha:g}_inner_record_consistent", flag(solve, "converged")) + raw, actual, correction, coefficient = (number(row, key) for key in ("raw_gravity_remainder", "corrected_gravity_remainder", "gravity_correction_norm", "correction_over_alpha_squared")) + floor = roundoff(*(number(row, key) for key in ("base_residual_norm", "raw_residual_norm", "corrected_residual_norm", "linear_prediction_norm"))) + self.check(f"curvature_{alpha:g}_finite_nonnegative", all(nonnegative(x) for x in (raw, actual, correction, coefficient))) + self.check(f"curvature_{alpha:g}_retained_unchanged", number(row, "retained_state_change") == 0.0) + self.check(f"curvature_{alpha:g}_coefficient_arithmetic", math.isclose(coefficient, correction / alpha**2, rel_tol=1e-12, abs_tol=1e-30)) + tolerance = number(solve, "threshold") + floor + self.check(f"curvature_{alpha:g}_remainder_cancelled", nonnegative(actual) and actual <= tolerance, + corrected_remainder=actual, allowed=tolerance, raw_remainder=raw) + if alpha <= .1 and raw > 100.0 * floor and correction > 100.0 * floor and number(solve, "threshold") < .01 * raw: + reliable.append((alpha, correction, coefficient)) + self.outcomes["curvature"] = {"fully_evaluated_steps": evaluated, "unavailable_steps": len(rows) - evaluated, + "resolved_small_steps": reliable} + if evaluated == 0: + self.incomplete("curvature_cancellation_available", reason="No prepared, solved, corrected trial") + reliable.sort(reverse=True) + if len(reliable) >= 2: + coarse, fine = reliable[-2:] + order = math.log(fine[1] / coarse[1]) / math.log(fine[0] / coarse[0]) + self.check("curvature_at_least_second_order_trend", order >= 1.7, observed_order=order, + minimum_test_order=1.7, alphas=[coarse[0], fine[0]], + interpretation="A finite-range trend check, not an asymptotic proof") + else: + self.incomplete("curvature_second_order_trend", reason="Fewer than two small steps above the heuristic arithmetic/solve floor") + + def reduced(self, summary, lookup): + tolerance = summary["gravity_completion_tolerance"] + self.check("gravity_completion_tolerance_valid", math.isfinite(tolerance) and tolerance > 0.0) + completions = read_rows(self.folder / "numerical_strategy_completions.csv", ("kind", "alpha", "completed", "final_gravity_norm", "retained_state_change")) + for index, row in enumerate(completions): + self.check(f"completion_{index}_retained_unchanged", number(row, "retained_state_change") == 0.0) + if flag(row, "completed"): + residual = number(row, "final_gravity_norm") + self.check(f"completion_{index}_prepared_gravity_gate", nonnegative(residual) and residual <= tolerance, + final_gravity_norm=residual, tolerance=tolerance) + self.outcomes["gravity_completion"] = {"attempts": len(completions), "completed": sum(flag(row, "completed") for row in completions)} + if summary.get("reduced_base_completed") != 1.0: + self.incomplete("reduced_base_available", reason="Base gravity completion did not reach its prepared-residual tolerance") + return + for kind in ("reduced_action", "reduced_lift"): + if (kind, 0.0, 0) not in lookup or not flag(lookup[kind, 0.0, 0], "converged"): + self.incomplete("reduced_action_available", reason=f"{kind} did not converge") + return + action = summary["reduced_action_norm"] + lifted_g = summary["lifted_gravity_action_norm"] + agreement = summary["lifted_retained_action_difference"] + lift = lookup["reduced_lift", 0.0, 0] + floor = roundoff(action, number(lift, "right_hand_side_norm")) + self.check("lifted_tangent_gravity_gate", nonnegative(lifted_g) and lifted_g <= number(lift, "threshold") + floor, + lifted_gravity_action=lifted_g, allowed=number(lift, "threshold") + floor) + allowed_agreement = 1e-8 * abs(action) + floor + self.check("schur_matches_lifted_full_action", nonnegative(agreement) and agreement <= allowed_agreement, + difference=agreement, allowed=allowed_agreement) + differences = read_rows(self.folder / "numerical_strategy_reduced_derivative.csv", ("alpha", "plus_completed", "minus_completed", "relative_error", "difference_quotient_error", "plus_gravity_norm", "minus_gravity_norm", "residual_subtraction_roundoff_scale")) + reliable = [] + valid = [] + for row in differences: + if not (flag(row, "plus_completed") and flag(row, "minus_completed")): + self.incomplete(f"reduced_fd_{row['alpha']}_completed_pair_available", + plus_completed=flag(row, "plus_completed"), minus_completed=flag(row, "minus_completed")) + continue + alpha, error, relative, floor = (number(row, key) for key in ("alpha", "difference_quotient_error", "relative_error", "residual_subtraction_roundoff_scale")) + valid.append((alpha, error, relative, floor)) + self.check(f"reduced_fd_{alpha:g}_finite", all(nonnegative(x) for x in (error, relative, floor))) + self.check(f"reduced_fd_{alpha:g}_completed_sides", all(nonnegative(number(row, key)) and number(row, key) <= tolerance for key in ("plus_gravity_norm", "minus_gravity_norm"))) + if error > 20.0 * floor: + reliable.append((alpha, error)) + self.outcomes["reduced_central_difference"] = {"completed_pairs": len(valid), "measurements": valid} + if not valid: + self.incomplete("reduced_central_difference_available", reason="No pair completed on both sides") + return + self.check("reduced_central_difference_accuracy", min(item[2] for item in valid) <= 1e-5, + minimum_relative_error=min(item[2] for item in valid), tolerance=1e-5) + reliable.sort(reverse=True) + orders = [math.log(b[1] / a[1]) / math.log(b[0] / a[0]) for a, b in zip(reliable, reliable[1:]) if a[0] != b[0]] + if any(order >= 1.5 for order in orders): + self.check("reduced_central_difference_second_order_trend", True, observed_orders=orders) + elif all(item[1] <= 20.0 * item[3] for item in valid): + self.check("reduced_central_difference_floor_consistent", True, + interpretation="Errors are at the recorded residual-subtraction scale; order is unresolved") + else: + self.incomplete("reduced_central_difference_trend", observed_orders=orders, + reason="No resolved second-order interval; inner completion error has no measured retained-residual amplification bound") + + +def assess(folder, required=(), expect_negative=False, restoration_tolerance=0.0): + assessment = Assessment(folder, restoration_tolerance) + sections = set(required) + if (folder / "polynomial_certificates.csv").exists() or expect_negative: + sections.add("polynomial") + if (folder / "numerical_strategy_summary.csv").exists(): + sections.add("gravity") + if not sections: + assessment.incomplete("probe_artifacts_available", reason="No recognized frozen strategy outputs") + for section in sorted(sections): + try: + if section == "polynomial": + assessment.polynomial(expect_negative) + else: + assessment.gravity() + except FileNotFoundError as error: + assessment.incomplete(f"{section}_artifacts_complete", error=str(error)) + except KeyError as error: + assessment.incomplete(f"{section}_records_complete", missing=str(error)) + except (ValueError, OSError, csv.Error, ZeroDivisionError) as error: + assessment.check(f"{section}_artifacts_valid", False, error=str(error)) + statuses = {item["status"] for item in assessment.checks} + status = "fail" if "fail" in statuses else "incomplete" if "incomplete" in statuses else "pass" + report = {"schema_version": 1, "result_folder": str(folder.resolve()), "status": status, + "checks": assessment.checks, "outcomes": assessment.outcomes, + "scope": "Frozen algebra, admissibility, and restoration checks; use check_plateau_physics.py for equilibrium/virial acceptance"} + return clean_json(report), {"pass": 0, "incomplete": 2, "fail": 3}[status] + + +def main(): + parser = argparse.ArgumentParser(description=__doc__) + parser.add_argument("result_folder", type=Path) + parser.add_argument("--require", choices=("polynomial", "gravity"), action="append", default=[]) + parser.add_argument("--expect-negative-base", action="store_true") + parser.add_argument("--restoration-tolerance", type=float, default=0.0) + args = parser.parse_args() + if not nonnegative(args.restoration_tolerance): + parser.error("--restoration-tolerance must be finite and nonnegative") + report, code = assess(args.result_folder, args.require, args.expect_negative_base, args.restoration_tolerance) + print(json.dumps(report, indent=2, allow_nan=False)) + return code + + +if __name__ == "__main__": + sys.exit(main()) diff --git a/experiments/check_numerical_strategies_README.md b/experiments/check_numerical_strategies_README.md new file mode 100644 index 0000000..9d00108 --- /dev/null +++ b/experiments/check_numerical_strategies_README.md @@ -0,0 +1,72 @@ +# Checking numerical-strategy experiments + +Run the checker on each completed frozen-probe directory: + +```sh +python3 experiments/check_numerical_strategies.py RESULT_FOLDER --require gravity +python3 experiments/check_numerical_strategies.py RESULT_FOLDER --require polynomial +python3 experiments/check_numerical_strategies.py LATE_ROTATING_FOLDER --require polynomial --expect-negative-base +python3 experiments/test_check_numerical_strategies.py +``` + +Both sections are discovered automatically when their CSV files exist. Repeating +`--require` demands both sections. Exit 0 means the requested checks passed; +3 means an invariant or numerical expectation failed; 2 means evidence is +incomplete. An honestly reported failed inner solve, unavailable trial, or +unresolved certificate must not turn into a passing verification. Expected +negative-base detection and rejection are successful guard tests, not successful +equilibrium solves. + +Accepted state, residual, physical state, saved correction, and retained-coordinate +changes must be exactly zero. Restored residual drift defaults to exactly zero; +`--restoration-tolerance` relaxes only that comparison when explicitly justified. +The JSON output records the chosen tolerance. + +For successful curvature probes, the checker compares the fresh gravity remainder +with the selected inner threshold plus a documented arithmetic scale, verifies +that only gravity changed, and checks correction/alpha² arithmetic. It requires +two small steps resolved above that scale to assess a second-order trend. An +observed order of at least 1.7 passes the finite-range trend check; this does not +prove an asymptotic rate. Higher order is allowed because a quadratic coefficient +can vanish. + +Reduced probes must satisfy the freshly prepared gravity completion tolerance, +the lifted gravity tangent gate, and agreement of the Schur action with the +retained part of the full lifted action. Central differences must reach relative +error at most 1e-5 and show a resolved interval with order at least 1.5, or reach +the recorded residual-subtraction scale. A plateau above that scale is incomplete +evidence rather than an automatic claim of a bad Jacobian: the retained-residual +amplification of the remaining gravity error has not been bounded. The recorded +roundoff scales are heuristics, not rigorous forward-error bounds. + +This script deliberately does not duplicate physical acceptance. For each +end-to-end model, run `check_plateau_physics.py` on its endpoint. That separate +gate requires recorded strict nonlinear convergence and all eight defined virial +measurements at or below 1e-5. Its local EOS, hydrostatic, and sign measurements +remain part of the physical assessment. + +## Bounded validation matrix + +These are proposed operational limits, not measured performance predictions. +Use one numerical process at a time and preserve all capped failures as results. + +| Layer | Cases | Bound / acceptance | +|---|---|---| +| Algebra and lifecycle | Gravity subsystem, Schur reduction, polynomial certificate tests; one and two MPI ranks | 2 minutes per test group; true-defect, stale-state, empty-local, rollback, type-contract coverage | +| Python checker | Synthetic passing, false-convergence, mutation, missing-data, negative-base, and NaN records | 30 seconds; standard-library unit tests | +| Positivity frozen states | n=1 J=0.1 captured iteration 2 and iteration 11; a certified seed/control | 10 minutes per probe; iteration 11 must be negative and reject repair with `invalid_base`; no claim that the earlier base is positive before checking it | +| Gravity curvature | n=1.5 J=0 captured iteration 1 | 15 minutes; alpha = 1, 0.1, 0.01, 0.001; gate fresh remainder cancellation and small-step order | +| Reduced gravity | Same n=1.5 state, plus constrained n=3 saved seed direction | 15 minutes per probe; completion and three central-difference pairs; n=3 direction remains an approximate diagnostic correction | +| Regression equilibrium | n=1 J=0 and low-rotation J=0.05 baseline; opt-in strategies individually | Preserve existing nonlinear/linear caps, add a 20-minute wall limit per model; strict solver and all virial gates remain mandatory | +| Failure-prone equilibrium | n=1 J=0.1 and n=1.5 J=0, strategies separately before combinations | 8 Newton steps / 20 minutes; report residual progress, auxiliary work, rejection causes, and physical endpoint without equating the cap with convergence | + +At the supplied frozen-gravity defaults, there are at most 20 explicitly requested +inner solves, each capped at 160 iterations, and at most 36 preparations in the +normal completed path. The maximum covers four curvature corrections, seven +gravity completions with up to two corrections each, and two reduced actions. +Rejected preparations can reduce those counts. Root integration controls and +end-to-end outer solves are separately bounded; their nested work is not covered +by this frozen-probe count. + +Treat outer reduced-solve time as inclusive of nested gravity solves. Report the +inner counters separately, but do not add the two times as if they were disjoint. diff --git a/experiments/test_check_numerical_strategies.py b/experiments/test_check_numerical_strategies.py new file mode 100644 index 0000000..5ada8bf --- /dev/null +++ b/experiments/test_check_numerical_strategies.py @@ -0,0 +1,170 @@ +#!/usr/bin/env python3 +"""Exercise checker acceptance and failure paths with synthetic probe records.""" + +import csv +import json +from pathlib import Path +import tempfile +import unittest + +from check_numerical_strategies import assess + + +def write_rows(folder, name, rows, fields=None): + fields = fields or list(rows[0]) + with (folder / name).open("w", newline="") as stream: + writer = csv.DictWriter(stream, fieldnames=fields) + writer.writeheader() + writer.writerows(rows) + + +def write_metrics(folder, name, values): + write_rows(folder, name, [{"metric": key, "value": value} for key, value in values.items()]) + + +def gravity_fixture(folder): + summary = {"accepted_state_change": 0, "accepted_residual_change": 0, "accepted_physical_change": 0, + "maximum_original_restoration_error": 0, "gravity_completion_tolerance": 1e-11, + "reduced_base_completed": 1, "reduced_action_norm": 2, + "lifted_gravity_action_norm": 1e-13, "lifted_retained_action_difference": 1e-13} + inner = [] + curvature = [] + for alpha in (1., .1, .01, .001): + rhs = .1 * alpha**2 + threshold = max(1e-13, 1e-9 * rhs) + residual = threshold * .1 + inner.append({"kind": "curvature", "alpha": alpha, "index": 0, "status": 0, "converged": 1, + "right_hand_side_norm": rhs, "true_residual_norm": residual, + "relative_true_residual_norm": residual / rhs, "threshold": threshold}) + curvature.append({"alpha": alpha, "prepared": 1, "inner_converged": 1, "corrected_prepared": 1, + "base_residual_norm": 1, "raw_residual_norm": 1, "corrected_residual_norm": .99, + "linear_prediction_norm": .99, "raw_gravity_remainder": rhs, + "corrected_gravity_remainder": residual, "gravity_correction_norm": .2 * alpha**2, + "correction_over_alpha_squared": .2, "retained_state_change": 0}) + for kind in ("reduced_action", "reduced_lift"): + inner.append({"kind": kind, "alpha": 0, "index": 0, "status": 0, "converged": 1, + "right_hand_side_norm": 1, "true_residual_norm": 1e-13, + "relative_true_residual_norm": 1e-13, "threshold": 1e-9}) + completions = [{"kind": "reduced_base", "alpha": 0, "completed": 1, + "final_gravity_norm": 5e-14, "retained_state_change": 0}] + differences = [] + for alpha in (.01, .001, .0001): + for kind in ("reduced_plus", "reduced_minus"): + completions.append({"kind": kind, "alpha": alpha, "completed": 1, + "final_gravity_norm": 5e-14, "retained_state_change": 0}) + differences.append({"alpha": alpha, "plus_completed": 1, "minus_completed": 1, + "relative_error": .25 * alpha**2, "difference_quotient_error": .5 * alpha**2, + "plus_gravity_norm": 5e-14, "minus_gravity_norm": 5e-14, + "residual_subtraction_roundoff_scale": 1e-14 / alpha}) + write_metrics(folder, "numerical_strategy_summary.csv", summary) + write_rows(folder, "numerical_strategy_inner.csv", inner) + write_rows(folder, "numerical_strategy_curvature.csv", curvature) + write_rows(folder, "numerical_strategy_completions.csv", completions) + write_rows(folder, "numerical_strategy_reduced_derivative.csv", differences) + return summary, inner, curvature, completions, differences + + +def polynomial_fixture(folder): + certificates = [{"state": name, "status": status, "witness_element": 459 if status == "negative" else -1, + "x": .5, "y": .5, "z": .5} + for name, status in (("frozen_base", "negative"), ("context_seed_anchor", "certified_nonnegative"), + ("original_full_candidate", "negative"))] + repair = [{"status": "invalid_base", "minimum_safe_fraction": 1, "base_status": "negative"}] + trials = [{"fraction": 0, "alpha": 1, "certificate": "negative", "prepared": 1, "merit_slope": -.1}] + write_rows(folder, "polynomial_certificates.csv", certificates) + write_rows(folder, "polynomial_repair.csv", repair) + write_rows(folder, "polynomial_repair_trials.csv", trials) + write_metrics(folder, "polynomial_restoration.csv", { + "accepted_physical_difference": 0, "accepted_normalized_difference": 0, "accepted_residual_difference": 0, + "restored_correction_difference": 0, "restored_residual_difference": 0}) + return certificates, repair, trials + + +class NumericalStrategyChecks(unittest.TestCase): + def setUp(self): + self.temporary = tempfile.TemporaryDirectory() + self.folder = Path(self.temporary.name) + + def tearDown(self): + self.temporary.cleanup() + + def test_complete_gravity_evidence_passes(self): + gravity_fixture(self.folder) + report, code = assess(self.folder, required=("gravity",)) + self.assertEqual(code, 0, report) + + def test_negative_late_base_is_an_expected_guard_success(self): + polynomial_fixture(self.folder) + report, code = assess(self.folder, expect_negative=True) + self.assertEqual(code, 0, report) + + def test_negative_base_cannot_be_repaired_silently(self): + _, repair, _ = polynomial_fixture(self.folder) + repair[0]["status"] = "feasible" + write_rows(self.folder, "polynomial_repair.csv", repair) + report, code = assess(self.folder, expect_negative=True) + self.assertEqual(code, 3, report) + + def test_false_convergence_report_fails(self): + _, inner, _, _, _ = gravity_fixture(self.folder) + inner[0]["true_residual_norm"] = 1.0 + write_rows(self.folder, "numerical_strategy_inner.csv", inner) + report, code = assess(self.folder) + self.assertEqual(code, 3, report) + + def test_changed_retained_state_fails(self): + _, _, curvature, _, _ = gravity_fixture(self.folder) + curvature[0]["retained_state_change"] = 1e-8 + write_rows(self.folder, "numerical_strategy_curvature.csv", curvature) + report, code = assess(self.folder) + self.assertEqual(code, 3, report) + + def test_failed_numerical_probe_is_incomplete_not_a_false_pass(self): + summary, inner, curvature, completions, _ = gravity_fixture(self.folder) + summary["reduced_base_completed"] = 0 + for row in inner: + row["status"], row["converged"] = 1, 0 + for row in curvature: + row["inner_converged"], row["corrected_prepared"] = 0, 0 + for row in completions: + row["completed"] = 0 + write_metrics(self.folder, "numerical_strategy_summary.csv", summary) + write_rows(self.folder, "numerical_strategy_inner.csv", inner) + write_rows(self.folder, "numerical_strategy_curvature.csv", curvature) + write_rows(self.folder, "numerical_strategy_completions.csv", completions) + report, code = assess(self.folder) + self.assertEqual(code, 2, report) + self.assertEqual(report["outcomes"]["gravity_inner"]["failed_calls"], 6) + + def test_missing_probe_is_incomplete(self): + report, code = assess(self.folder, required=("gravity",)) + self.assertEqual(code, 2, report) + + def test_unfinished_restoration_record_is_incomplete(self): + summary, _, _, _, _ = gravity_fixture(self.folder) + del summary["accepted_state_change"] + write_metrics(self.folder, "numerical_strategy_summary.csv", summary) + report, code = assess(self.folder) + self.assertEqual(code, 2, report) + + def test_nonfinite_restoration_fails_and_serializes_standard_json(self): + summary, _, _, _, _ = gravity_fixture(self.folder) + summary["accepted_state_change"] = float("nan") + write_metrics(self.folder, "numerical_strategy_summary.csv", summary) + report, code = assess(self.folder) + self.assertEqual(code, 3, report) + json.dumps(report, allow_nan=False) + + def test_only_restoration_not_accepted_state_may_use_optional_tolerance(self): + summary, _, _, _, _ = gravity_fixture(self.folder) + summary["maximum_original_restoration_error"] = 1e-13 + write_metrics(self.folder, "numerical_strategy_summary.csv", summary) + self.assertEqual(assess(self.folder)[1], 3) + self.assertEqual(assess(self.folder, restoration_tolerance=1e-12)[1], 0) + summary["accepted_state_change"] = 1e-13 + write_metrics(self.folder, "numerical_strategy_summary.csv", summary) + self.assertEqual(assess(self.folder, restoration_tolerance=1e-12)[1], 3) + + +if __name__ == "__main__": + unittest.main()