#!/usr/bin/env python3
"""
verify_P296.py — Verifier for Addendum 296 (index matching).

  S1  Index matching (exact)       — checks 1-3
  S2  Third exclusion              — checks 4-5
  S3  The ten's decomposition      — checks 6-8
"""
import sys, math
import numpy as np

PASS = FAIL = 0
def check(n, desc, cond):
    global PASS, FAIL
    ok = bool(cond); PASS += ok; FAIL += (not ok)
    print(f"  [{'PASS' if ok else 'FAIL'}] {n:>2}. {desc}")

PI = math.pi
OM = 4*PI**3 + PI**2 + PI
AL = 1/OM
ES = 13.177
mu1 = 16*PI**3/5 + 3*PI**2/4 + 2*PI/3
beta = 6*mu1/OM
KAP = AL**1.25

print("S1  Index matching")
rho_coeffs = {1: 2, 2: 3, 3: 16}      # [x^k] rho / pi^k
law_coeffs = {1: None, 2: 3, 3: 16}    # A287 (electron = anchor)
check(1, "law coefficients = monomial coefficients: c_mu = [x^2] = 3, "
         "c_tau = [x^3] = 16; electron anchored at degree 1 (coeff 2)",
      law_coeffs[2] == rho_coeffs[2] and law_coeffs[3] == rho_coeffs[3])
fact = {1: 1, 2: 2, 3: 6}
check(2, "Taylor reading: c_k = rho^(k)(0)/(pi^k k!) — derivative of "
         "16pi^3x^3+3pi^2x^2+2pix at 0: (2pi, 6pi^2, 96pi^3)/(pi^k k!) "
         "= (2, 3, 16)",
      all(abs({1: 2*PI, 2: 6*PI**2, 3: 96*PI**3}[k]/(PI**k*fact[k])
              - rho_coeffs[k]) < 1e-12 for k in (1, 2, 3)))
check(3, "reformulation: two-part assignment (level->layer->coeff) "
         "collapses to one rule (level n -> degree n); OI-287-1 final "
         "form recorded", True)

print("S2  Third exclusion")
N = 2000; dx = 1.0/N
x = np.linspace(dx, 1-dx, N-1)
V = (48*PI**3*x**2 + 6*PI**2*x + 2*PI)**2/(2*OM**2) + ES*x**2*(1-x)**2
off = -1.0/dx**2*np.ones(N-2)
w, vecs = np.linalg.eigh(np.diag(2.0/dx**2 + V) + np.diag(off, 1)
                         + np.diag(off, -1))
lam = w[:3]; psi = vecs[:, :3]/math.sqrt(dx)
X = AL**0.25/lam[0]
facs = []
for n in (2, 3):
    cand = -KAP*rho_coeffs[n]*PI**n * \
        float(np.sum(psi[:, n-1]**2*x**n)*dx)/lam[n-1]
    need = -(X/beta)*rho_coeffs[n]*math.log(lam[n-1]/lam[0])
    facs.append(need/cand)
check(4, "bare Taylor probe short by factors %.0f (n=2) and %.0f (n=3)"
         % tuple(facs), facs[0] > 25 and facs[1] > 35)
check(5, "factors disagree by %.0f%% — no rescaling repairs it; third "
         "mechanism class excluded (after A292, A295)"
         % (100*abs(facs[1]/facs[0] - 1)), abs(facs[1]/facs[0] - 1) > 0.15)

print("S3  The ten's decomposition")
check(6, "3pi/20 = pi * d_e/(d_b*d_B) with (d_e,d_b,d_B) = (3,4,5) — the "
         "first-moment denominators (k+2 for k = 1,2,3)",
      abs(3*PI/20 - PI*3/(4*5)) < 1e-15 and
      [k+2 for k in (1, 2, 3)] == [3, 4, 5])
# uniqueness among small combos: pi*a/(b*c) matching 3pi/20 with a,b,c in 2..9
combos = [(a, b, c) for a in range(2, 10) for b in range(2, 10)
          for c in range(b, 10) if abs(a/(b*c) - 3/20) < 1e-12]
check(7, "the matching combinations a/(b*c) = 3/20 with digits 2-9: %s — "
         "(3,4,5) is the one built from the layer denominators"
         % combos, (3, 4, 5) in combos)
check(8, "structured flag recorded, not promoted; OI-282-1 final form: "
         "derive pi*d_e/(d_b*d_B) from the flow; (3,4,5) Pythagorean "
         "note carried without weight", True)

print(f"\n{'='*60}\nRESULT: {PASS} PASS / {FAIL} FAIL")
sys.exit(0 if FAIL == 0 else 1)