#!/usr/bin/env python3
"""
verify_P080.py -- Addendum 80: MU from the spectral mean of rho.

This verifier checks the rho-moment arithmetic and the claimed closure of the
P77 spectral-angle assumption.  The moment identity itself is correct:
MU=mu1/mu0=CGF'(0) is the rho-weighted mean of x.  The flagged issues are the
stale printed MU value and proof/status jumps: the Hopf derivative calculation
does not equal the rho mean, the mass-coordinate conversion direction is
reversed in one sentence, and the paper's unconditional Phase-5b closure rests
on that unproved angular-conversion identification.
"""

from __future__ import annotations

import math
import sys
from pathlib import Path

import mpmath as mp

sys.path.insert(0, str(Path(__file__).resolve().parent))
from verify_common import CheckResult, Verifier


class ModernVerifier(Verifier):
    """Output adapter: modern check-line format. Tolerance logic is
    inherited unchanged from verify_common.Verifier; only printing and
    the footer differ."""

    def __init__(self, name: str):
        self.name = name
        self.results = []
        print(name)

    def record(self, label, ok, computed="", claimed="", detail=""):
        self.results.append(CheckResult(label, ok, computed, claimed, detail))
        n = len(self.results)
        expected = (not ok) and ("Expected" in str(detail))
        desc = f"{label} -- {detail}" if expected else label
        print(f"  [{'PASS' if ok else 'FAIL'}] {n:>2}. {desc}")
        if computed != "" or claimed != "":
            print(f"       computed: {computed}")
            print(f"       claimed : {claimed}")
        if detail and not expected:
            print(f"       {detail}")
        return ok

    def summary(self):
        passed = sum(r.ok for r in self.results)
        failed = len(self.results) - passed
        print(f"\n{'='*60}\nRESULT: {passed} PASS / {failed} FAIL")
        return 1 if failed else 0



mp.mp.dps = 60
v = ModernVerifier("P080 -- MU as Spectral Mean of rho")

ROOT = Path(__file__).resolve().parents[1]
TEX = (ROOT / "80_Addendum_MUFromDensity.tex").read_text()

PI = mp.pi
MU0 = 4 * PI**3 + PI**2 + PI
MU1 = 16 * PI**3 / 5 + 3 * PI**2 / 4 + 2 * PI / 3
MU = MU1 / MU0
EPS = MU * PI / 6
EPL = mp.mpf("68.096")
ER_PRED = EPL - 12 + EPS
ER_OBS = mp.mpf("56.502")
ER_SIGMA = mp.mpf("0.009")
EV = mp.mpf("19.63556749202108")
ME_EV = mp.mpf("511000")


def f(x: mp.mpf) -> float:
    return float(x)


def rho(x: mp.mpf) -> mp.mpf:
    return 2 * PI * x + 3 * PI**2 * x**2 + 16 * PI**3 * x**3


v.check("mu0 exact integral", f(MU0), 137.036, rel=3e-6)
v.check("mu1 exact integral", f(MU1), 108.717, rel=5e-6)
v.check("mu1 edge component", f(2 * PI / 3), 2.0944, rel=2e-5)
v.check("mu1 boundary component", f(3 * PI**2 / 4), 7.4022, rel=2e-5)
v.check("mu1 bulk component", f(16 * PI**3 / 5), 99.220, rel=5e-5)
v.check(
    "printed MU numeric",
    f(MU),
    0.793338,
    rel=1e-6,
    detail="Expected fail: exact mu1/mu0≈0.793342208, so 0.793338 is stale at six decimals.",
)
v.check("MU rounded corpus value", f(MU), 0.79334, rel=5e-6)
v.check("MU reduced formula", f((192 * PI**2 + 45 * PI + 40) / (240 * PI**2 + 60 * PI + 60)), f(MU), rel=1e-12)
v.check("CGF'(0)=mu1/mu0", f(MU1 / MU0), f(MU), rel=1e-12)
v.check("spectral mean lies below bulk 4/5", f(MU), 0.8, rel=9e-3)

v.check("sector ratio EB/mu0", f((4 * PI**3) / MU0), 0.9051, rel=7e-5)
v.check("sector ratio (Eb+EB)/mu0", f((PI**2 + 4 * PI**3) / MU0), 0.9771, rel=3e-5)
v.check("alpha-minus-one ratio", f((MU0 - 1) / MU0), 0.9927, rel=5e-5)
v.check("boundary EW ratio", f(PI**2 / (PI + PI**2)), 0.7585, rel=7e-5)
v.check("mu1 sector decomposition", f((mp.mpf(2) / 3) * PI + (mp.mpf(3) / 4) * PI**2 + (mp.mpf(4) / 5) * 4 * PI**3), f(MU1), rel=1e-12)

v.check("epsilon MU*pi/6", f(EPS), 0.415392, rel=5e-6)
v.check("E_Pl - E_R predicted gap", f(12 - EPS), 11.584608, rel=5e-6)
v.check("E_R predicted", f(ER_PRED), 56.511, rel=8e-6)
v.check("E_R sigma with displayed uncertainty", f(abs(ER_PRED - ER_OBS) / ER_SIGMA), 1.02, rel=3e-2)

v.record(
    "mass-coordinate conversion sentence uses the right direction",
    False,
    computed="from E=pi+ln(m/me)/MU, dE/dln(m/me)=1/MU",
    claimed="MU is the conversion factor from log-mass to spectral units",
    detail="Expected audit fail: the written derivative gives 1/MU for log-mass -> E; MU is the inverse conversion.",
)

avg_abs_derivative = mp.quad(lambda x: 2 * mp.sqrt(x * (1 - x)) * rho(x), [0, 1]) / MU0
derivative_at_mean = 2 * mp.sqrt(MU * (1 - MU))
v.check(
    "Hopf average |dx/dbeta| equals MU",
    f(avg_abs_derivative),
    f(MU),
    rel=5e-3,
    detail="Expected fail: averaging the actual |d cos^2(beta)/dbeta| gives about 0.693, not MU.",
)
v.check(
    "Hopf derivative at mean x equals MU",
    f(derivative_at_mean),
    f(MU),
    rel=5e-3,
    detail="Expected fail: evaluating the derivative magnitude at x=MU gives about 0.810, not MU.",
)
v.record(
    "Assumption 4.1 follows from the moment computation alone",
    False,
    computed="moment computation proves MU=<x>_rho",
    claimed="one radian of G2 Cartan angle corresponds to MU spectral units",
    detail="Expected proof-audit fail: the angular-conversion map is identified with the mean, but not derived from the Hopf/G2 maps.",
)
v.record(
    "Phase 5b is closed without residual assumptions",
    False,
    computed="closure depends on the unproved angular-conversion identification above",
    claimed="Phase 5b closed without qualification",
    detail="Expected status fail: the rho mean is proved, but the conversion-factor theorem is not established by the displayed calculation.",
)

enu = 2 * EV - ER_PRED
mnu_mev = ME_EV * mp.e ** (MU * (enu - PI)) * 1000
v.check("neutrino energy from exact rounded inputs", f(enu), -17.239, rel=1e-4)
v.check(
    "predicted neutrino mass",
    f(mnu_mev),
    48.7,
    rel=2e-3,
    detail="Expected fail: exact use of the displayed inputs gives about 48.53 meV.",
)
v.check(
    "neutrino sigma vs 49.5±0.35 meV",
    f(abs(mp.mpf("49.5") - mnu_mev) / mp.mpf("0.35")),
    2.3,
    rel=8e-2,
    detail="Expected fail: 49.5±0.35 against 48.53 meV is about 2.8 sigma.",
)
v.record(
    "neutrino prediction uses no oscillation input",
    "No oscillation data enters" in TEX,
    computed="explicit no-oscillation-data claim present",
    claimed="status wording recorded",
)

sys.exit(v.summary())
