diff --git a/experiments/audit_density_order.py b/experiments/audit_density_order.py index dac3d22..d208eb9 100644 --- a/experiments/audit_density_order.py +++ b/experiments/audit_density_order.py @@ -609,11 +609,24 @@ def audit_case(campaign, name): "linear_iterations": sum(int(r["linear_iterations"]) for r in history)} fine = totals[("stored_density", degree, 30)] physical18 = result["physical"]["18"]["metrics"] - for region_key, physical_key in (("mass", "mass"), ("volume", "volume"), ("inertia_z", "moment_of_inertia_z")): + for region_key, physical_key in (("mass", "mass"), ("volume", "volume")): close(fine[region_key], physical18[physical_key], "Sliced/physical " + region_key, rtol=2e-9, atol=1e-13) + # Position factors raise the normal polynomial degree beyond the + # eight-point sliced rule's exactness. Tangential refinement does not + # test that error, so moment differences are observations, not identities. + radius = float(meta["R_reference"]) + require(math.isfinite(radius) and radius > 0, "Invalid reference radius") + moment_comparisons = {"scope": "Order14 normal slicing versus unsliced q18; first moments and inertia may have normal quadrature error."} + values = [("inertia_z", physical18["moment_of_inertia_z"], 2)] for j, axis in enumerate("xyz"): - close(fine[f"first_moment_{axis}"], physical18[f"center_of_mass_{j}"] * physical18["mass"], - "Sliced/physical first moment", atol=2e-12) + values.append((f"first_moment_{axis}", physical18[f"center_of_mass_{j}"] * physical18["mass"], 1)) + for key, reference, power in values: + require(math.isfinite(fine[key]) and math.isfinite(reference), "Nonfinite moment observation") + difference = fine[key] - reference + mass_radius_scale = abs(physical18["mass"]) * radius ** power + moment_comparisons[key] = {"sliced": fine[key], "physical_q18": reference, "difference": difference, + "difference_over_mass_radius_power": difference / mass_radius_scale if mass_radius_scale else None} + result["sliced_vs_physical_moments"] = moment_comparisons return result @@ -733,6 +746,7 @@ def main(): "scope": __doc__, "limitations": ["Local projection M,b and mapped weights are not exported.", "Complete negative partitions are checked only on lines with exported negative intervals; no global positivity claim.", "Shallow unresolved witnesses remain recorded; root position tests have coefficient-scaled allowances.", + "Tangential convergence does not establish accuracy of the fixed eight-point normal rule; sliced first moments and inertia are compared observationally with q18.", "Coupled timing and convergence gates are producer observations, not independent equation evaluations."], "campaigns": {}, "consumed_sha256": {str(Path(common.__file__).resolve()): digest(Path(common.__file__))}} comparable = defaultdict(list)