#!/usr/bin/env python3
"""
verify_P007.py -- Paper 07: Self-Referential Observation / lepton masses.

This verifier checks the moment arithmetic, the observation-operator eigenvalue
numerics, and the power-law lepton mass predictions in
toe/07_Paper_SelfReferentialObservation.tex.  The broad mass-ratio numerics
reproduce at the paper's rough precision.  The flagged issues are that the
finite-difference matrix matches Dirichlet-like boundary conditions rather than
the stated Neumann condition at x=1, the Planck-log comparison is arithmetically
wrong, and the "zero free parameters"/mixing-angle claims are not supported by
the formulas printed in the paper.
"""

from __future__ import annotations

import math
import sys
from pathlib import Path

import numpy as np

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


class ModernVerifier(Verifier):
    """Local output adapter: inherits Verifier's tolerance logic unchanged,
    emits the corpus's modern check-line format (numbered [PASS]/[FAIL]
    lines, computed/claimed as indented info lines)."""

    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)
        desc = label if (ok or not detail) else f"{label} -- {detail}"
        print(f"  [{'PASS' if ok else 'FAIL'}] {n:>2}. {desc}")
        if computed != "" or claimed != "":
            print(f"        computed: {computed}")
            print(f"        claimed : {claimed}")
        if ok and detail:
            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


v = ModernVerifier("P007 -- Self-Referential Observation")

ROOT = Path(__file__).resolve().parents[2]
TEX = (ROOT / "toe" / "07_Paper_SelfReferentialObservation.tex").read_text()

PI = math.pi
MU0 = 4 * PI**3 + PI**2 + PI
MU1 = 16 * PI**3 / 5 + 3 * PI**2 / 4 + 2 * PI / 3
MU2 = 8 * PI**3 / 3 + 3 * PI**2 / 5 + PI / 2
E_SELF = 13.177
ME = 0.511
MMU_EXP = 105.7
MTAU_EXP = 1777.0


def pct(value: float, target: float) -> float:
    return 100 * (value - target) / abs(target)


def observation_eigenvalues(neumann_at_one: bool = False) -> np.ndarray:
    n = 1000
    dx = 1 / n
    x = np.linspace(dx, 1, n)
    rho_prime = 48 * PI**3 * x**2 + 6 * PI**2 * x + 2 * PI
    potential = rho_prime**2 / (2 * MU0**2) + E_SELF * x**2 * (1 - x) ** 2
    diag = np.full(n, 2 / dx**2) + potential
    if neumann_at_one:
        diag[-1] = 1 / dx**2 + potential[-1]
    off = np.full(n - 1, -1 / dx**2)
    matrix = np.diag(diag) + np.diag(off, 1) + np.diag(off, -1)
    return np.linalg.eigvalsh(matrix)[:4]


eigs_dirichlet_like = observation_eigenvalues(neumann_at_one=False)
eigs_mixed = observation_eigenvalues(neumann_at_one=True)

beta1 = 6 * MU1 / MU0
beta2 = beta1 * MU2 / MU1
lambda1, lambda2, lambda3 = 16.52, 51.28, 102.07
mu_pred = ME * (lambda2 / lambda1) ** beta1
tau_pred = mu_pred * (lambda3 / lambda2) ** beta2
avg_abs_error = (abs(pct(mu_pred, MMU_EXP)) + abs(pct(tau_pred, MTAU_EXP))) / 2


v.check("mu0", MU0, 137.036, rel=3e-6)
v.check("mu1", MU1, 108.717, rel=5e-6)
v.check("mu2", MU2, 90.176, rel=5e-7)
v.check("mu1/mu0", MU1 / MU0, 0.7933, rel=6e-5)
v.check("mu2/mu1", MU2 / MU1, 0.8295, rel=6e-5)
v.check("kappa=alpha^(5/4)", (1 / MU0) ** (5 / 4), 0.002133, rel=2e-3)
v.check("n0=80+pi", 80 + PI, 83.14, rel=2e-5)

v.check("beta empirical log(206.85)/log(3.10)", math.log(206.85) / math.log(3.10), 4.708, rel=2e-3)
v.check("5*mu1/mu0", 5 * MU1 / MU0, 3.967, rel=8e-5)
v.check("6*mu1/mu0", beta1, 4.760, rel=2e-5)
v.check("7*mu1/mu0", 7 * MU1 / MU0, 5.553, rel=8e-5)
v.check("4/(mu2/mu1)", 4 / (MU2 / MU1), 4.822, rel=9e-5)
v.check("beta2 recursive", beta2, 3.948, rel=8e-5)

v.check("lambda1 from printed matrix discretization", float(eigs_dirichlet_like[0]), 16.52, rel=2e-3)
v.check("lambda2 from printed matrix discretization", float(eigs_dirichlet_like[1]), 51.28, rel=2e-3)
v.check("lambda3 from printed matrix discretization", float(eigs_dirichlet_like[2]), 102.07, rel=2e-3)
v.record(
    "printed matrix implements Neumann boundary at x=1",
    False,
    computed=f"Dirichlet-like first eigenvalues {eigs_dirichlet_like[:3]}; mixed Dirichlet/Neumann values {eigs_mixed[:3]}",
    claimed="psi(0)=0 and dpsi/dx|_{x=1}=0",
    detail="Expected fail: the listed eigenvalues match the unmodified tridiagonal matrix, not a Neumann endpoint implementation.",
)
v.record(
    "x^2(1-x)^2 potential enforces the wavefunction boundary conditions",
    False,
    computed="the factor makes the potential term vanish at the endpoints",
    claimed="it enforces psi(0)=0 and psi'(1)=0",
    detail="Expected proof-audit fail: vanishing potential does not impose boundary conditions on psi.",
)

v.check("muon mass prediction", mu_pred, 112.0, rel=3e-3)
v.check("tau mass prediction", tau_pred, 1697.0, rel=2e-3)
v.check("muon/electron ratio", mu_pred / ME, 219.0, rel=3e-3)
v.check("tau/muon ratio", tau_pred / mu_pred, 15.1, rel=4e-3)
v.check("muon error percent", abs(pct(mu_pred, MMU_EXP)), 6.0, rel=3e-2)
v.check("tau error percent", abs(pct(tau_pred, MTAU_EXP)), 4.5, rel=4e-2)
v.check("quality score", 100 - avg_abs_error, 94.8, rel=1e-3)

planck_mev = 1.2209e22
planck_log_relation = 1 / math.log(planck_mev / ME)
v.check(
    "mu1/mu0 Planck-log comparison",
    planck_log_relation,
    MU1 / MU0,
    rel=1e-2,
    detail="Expected fail: the printed log(e)/log(M_Pl/m_e)/log(e) is about 0.0194, not 0.7933.",
)
v.record(
    "beta has no scan/selection input",
    False,
    computed="paper tests n=5,6,7 and selects n=6 because it is closest to the empirical beta",
    claimed="purely geometric with no fitting required",
    detail="Expected status fail: this is a small discrete scan unless the factor 6 is independently derived.",
)
v.record(
    "zero free parameters",
    False,
    computed="electron mass is used as an input; E_self is imported; beta factor is selected from a scan",
    claimed="zero free parameters",
    detail="Expected status fail.",
)

Y = np.array([[12.35, 0.87, 0.15], [0.87, 8.92, 1.45], [0.15, 1.45, 6.78]])
pair_angles = []
for i, j in [(0, 1), (0, 2), (1, 2)]:
    pair_angles.append(0.5 * math.degrees(math.atan2(2 * Y[i, j], Y[i, i] - Y[j, j])))
v.record(
    "printed Yukawa mixing angles follow from the displayed matrix",
    False,
    computed=f"simple Jacobi pair angles are {[round(a, 2) for a in pair_angles]} degrees",
    claimed="theta12≈62 deg, theta13≈31 deg, theta23≈84 deg",
    detail="Expected verification fail: no diagonalization convention is given that reproduces the printed angles from the displayed matrix.",
)
v.record(
    "three-family count follows from S3 topology as a theorem",
    False,
    computed="S^3 has three intrinsic dimensions, but no spectral/counting theorem is given tying this to exactly three fermion families",
    claimed="three families because S^3 has three dimensions",
    detail="Expected proof-audit fail.",
)
v.record(
    "Yukawa matrix entries are derivable from printed eigenfunctions",
    False,
    computed="the matrix is quoted without the eigenfunctions, quadrature code, or normalization convention needed to recompute it",
    claimed="computed from y_ij integral",
    detail="Expected verification gap.",
)

sys.exit(v.summary())