#!/usr/bin/env python3
# probe_selection_framing.py
# Copyright Léon Fernando Vlegels -- CC BY 4.0
#
# A321 framing audit for OI-287-1.
#
# Question: since A296/A299 the thread (A296-A320, ~18 probes) has
# operationalized the selection rule as an OVERLAP-DIAGONAL functional
#       B[n] = <psi_n | c_k pi^k x^k | psi_n> / lambda_n
# and asked whether B reproduces the target (c1,c2,c3)=(2,3,16),
# normalized (1, 1.5, 8). Seven arenas agree on overshoot-or-undershoot
# with nothing reaching mu=1.5. This script re-grounds the framing
# against the corpus text:
#
#   1. ORIGIN of (2,3,16): verify they are the Taylor coefficients of
#      rho(x) = 16 pi^3 x^3 + 3 pi^2 x^2 + 2 pi x, i.e. c_k = rho^(k)(0)/(pi^k k!).
#      The TARGET is literally the density's own defining coefficients.
#   2. THE CLAIM the corpus actually states (A287, A296): the family<->
#      coefficient relationship is a LEPTON MASS correction law
#      beta_n = 6 mu1/mu0 - c_n alpha^(1/4)/lambda_1, NOT an overlap-diagonal.
#   3. DRIFT CHECK: the overlap-diagonal <psi_n|c_k pi^k x^k|psi_n> is
#      A296's "bare Taylor probe", which A296 prop. (Third exclusion)
#      and A297 BOTH EXCLUDED. It was revived in A299 as a "lead" and
#      carried forward for ~18 probes.
#   4. ALTERNATIVE: the corpus-faithful object is the readout of the
#      coefficients themselves (Taylor data), which gives (2,3,16) /
#      (1,1.5,8) TRIVIALLY by definition; the magnitude "gap" is an
#      artifact of matching the overlap-diagonal, a functional the
#      corpus excluded. The corpus's OTHER use of rho is the moment
#      sequence mu_n = int x^n rho dx, used for the mass hierarchy
#      (P20 line 120,133-142); these are a monotone-decreasing
#      sequence and do NOT equal (2,3,16) either -- the target is the
#      coefficients, not the moments.

import json
import os
import numpy as np

PI = np.pi


def rho(x):
    return 16 * PI**3 * x**3 + 3 * PI**2 * x**2 + 2 * PI * x


def rho_prime(x):
    return 48 * PI**3 * x**2 + 6 * PI**2 * x + 2 * PI


def taylor_coeff(k):
    """c_k = rho^(k)(0) / (pi^k k!), in pi^k units."""
    # rho = 2 pi x + 3 pi^2 x^2 + 16 pi^3 x^3
    # rho'(0)   = 2 pi              -> c1 = 2 pi / (pi^1 1!)   = 2
    # rho''(0)  = 2 * 3 pi^2        -> c2 = 6 pi^2 / (pi^2 2!) = 3
    # rho'''(0) = 6 * 16 pi^3       -> c3 = 96 pi^3 / (pi^3 3!) = 16
    from math import factorial
    derivs = {1: 2 * PI, 2: 6 * PI**2, 3: 96 * PI**3}
    return derivs[k] / (PI**k * factorial(k))


def moment(n):
    """mu_n = int_0^1 x^n rho(x) dx."""
    return 2 * PI / (n + 2) + 3 * PI**2 / (n + 3) + 16 * PI**3 / (n + 4)


def build_H(N=1500):
    """Production operator fixed across all probes (OPEN_HEARTS_STATUS, A287)."""
    Omega = 4 * PI**3 + PI**2 + PI
    E_self = 13.177
    x = np.linspace(0.0, 1.0, N + 2)[1:-1]  # interior, Dirichlet
    h = 1.0 / (N + 1)
    V = (rho_prime(x))**2 / (2 * Omega**2) + E_self * x**2 * (1 - x)**2
    main = 2.0 / h**2 + V
    off = -1.0 / h**2 * np.ones(N - 1)
    H = np.diag(main) + np.diag(off, 1) + np.diag(off, -1)
    w, v = np.linalg.eigh(H)
    # normalize eigenvectors in L2 with trapezoid weight h
    for j in range(v.shape[1]):
        nrm = np.sqrt(np.sum(v[:, j]**2) * h)
        v[:, j] /= nrm
    return x, h, w[:3], v[:, :3]


def overlap_diagonal(x, h, lam, psi, c):
    """B[n] = <psi_n | c_n pi^n x^n | psi_n> / lambda_n  (A299 base)."""
    B = np.zeros(3)
    for n in range(3):
        k = n + 1
        op = c[n] * PI**k * x**k
        B[n] = np.sum(psi[:, n]**2 * op) * h / lam[n]
    return B


def main():
    out = {}

    # --- 1. ORIGIN: (2,3,16) are rho's Taylor coefficients ---
    c = [taylor_coeff(1), taylor_coeff(2), taylor_coeff(3)]
    out["taylor_coefficients"] = {"c1": c[0], "c2": c[1], "c3": c[2]}
    out["taylor_match_target"] = bool(
        abs(c[0] - 2) < 1e-12 and abs(c[1] - 3) < 1e-12 and abs(c[2] - 16) < 1e-12
    )
    # rho' integrates to rho with coefficients (2,3,16) at degrees 1,2,3
    out["rho_prime_to_rho_arithmetic"] = {
        "rho_prime": "48 pi^3 x^2 + 6 pi^2 x + 2 pi",
        "rho": "16 pi^3 x^3 + 3 pi^2 x^2 + 2 pi x",
        "coeffs_pi_units_deg_1_2_3": [2, 3, 16],
    }

    # TARGET readout: the coefficients themselves, normalized by c1.
    target = np.array([2.0, 3.0, 16.0])
    target_norm = target / target[0]
    out["target"] = target.tolist()
    out["target_normalized"] = target_norm.tolist()  # (1, 1.5, 8)
    # The corpus-faithful "functional" for "family n draws coefficient n"
    # is the direct coefficient readout: it equals the target trivially.
    coeff_readout = np.array(c)
    out["coefficient_readout"] = coeff_readout.tolist()
    out["coefficient_readout_normalized"] = (coeff_readout / coeff_readout[0]).tolist()
    out["coefficient_readout_L2_to_target"] = float(
        np.linalg.norm(coeff_readout / coeff_readout[0] - target_norm)
    )

    # --- 2/4. The moments mu_n: the corpus's OTHER use of rho (P20) ---
    mu = [moment(n) for n in range(4)]
    out["moments_mu"] = {f"mu_{n}": mu[n] for n in range(4)}
    out["mu0_equals_alpha_inv"] = bool(abs(mu[0] - (4 * PI**3 + PI**2 + PI)) < 1e-9)
    out["moment_ratios_over_mu0"] = [mu[n] / mu[0] for n in range(4)]
    out["moment_ratios_over_mu1"] = [mu[n] / mu[1] for n in range(4)]
    # do moments (or their ratios) reproduce (2,3,16)/(1,1.5,8)?
    mom_ratio_vec = np.array([mu[1] / mu[0], mu[2] / mu[0], mu[3] / mu[0]])
    out["moment_ratio_vec_n1_2_3"] = mom_ratio_vec.tolist()
    out["moments_give_target"] = bool(
        np.linalg.norm(mom_ratio_vec / mom_ratio_vec[0] - target_norm) < 0.1
    )
    out["moments_L2_to_target"] = float(
        np.linalg.norm(mom_ratio_vec / mom_ratio_vec[0] - target_norm)
    )

    # --- 3. DRIFT CHECK: the overlap-diagonal (A296 excluded, A299 revived) ---
    x, h, lam, psi = build_H()
    out["eigenvalues"] = lam.tolist()
    B = overlap_diagonal(x, h, lam, psi, c)
    Bn = B / B[0]
    out["overlap_diagonal_B"] = B.tolist()
    out["overlap_diagonal_normalized"] = Bn.tolist()
    out["overlap_diagonal_L2_to_target"] = float(np.linalg.norm(Bn - target_norm))
    out["overlap_diagonal_mu_entry"] = float(Bn[1])  # the stubborn ~1.12 entry

    # --- VERDICT ---
    # FRAMING-DRIFT: the overlap-diagonal is NOT mandated by the corpus.
    # A287/A296 state the rule as a lepton MASS correction law whose c_n
    # are the density's Taylor coefficients, read directly. The direct
    # coefficient readout gives (2,3,16)/(1,1.5,8) trivially (L2=0).
    # The overlap-diagonal is A296's "bare Taylor probe", EXCLUDED in
    # A296 and A297, revived as a "lead" in A299 and carried forward
    # ~18 probes. The magnitude "gap" (mu=1.5 unreached) is an artifact
    # of matching a functional the corpus excluded.
    drift = (
        out["taylor_match_target"]
        and out["coefficient_readout_L2_to_target"] < 1e-9
        and out["overlap_diagonal_L2_to_target"] > 0.1
    )
    out["verdict"] = "FRAMING-DRIFT" if drift else "INCONCLUSIVE"

    out["corpus_citations"] = {
        "origin_2_3_16": "287_Addendum_LeptonCorrectionLaw.tex lines 40-54, 105-119 "
        "(rho(x)=16 pi^3 x^3 + 3 pi^2 x^2 + 2 pi x; c_n the density coefficients)",
        "rule_is_mass_law": "287 lines 49-58 (beta_n = 6 mu1/mu0 - c_n alpha^(1/4)/lambda_1; "
        "reproduces m_mu/m_e, m_tau/m_e at 1e-4)",
        "taylor_reformulation": "296_Addendum_IndexMatching.tex lines 82-95 "
        "(c_k = rho^(k)(0)/(pi^k k!); selection rule = level n draws degree n)",
        "overlap_diagonal_excluded": "296 lines 97-107 (Third exclusion: bare Taylor probe "
        "diagonal coupling fails by factors 32,42) and 297 lines 61-77 (no index matching "
        "in position moments)",
        "overlap_diagonal_revived_as_lead": "299_Addendum_GapWeightLead.tex lines 46-55 "
        "(density-Taylor overlap / lambda_n carried forward as the 'strongest lead')",
        "moment_use_in_corpus": "20_Paper_StandardPhysicsEmbedding.tex lines 120,133-142 "
        "(mu_n = int x^n rho dx; mass hierarchies from mu_n; mu_0 = alpha^-1) and "
        "lines 451,459-461 (m_n ~ exp((mu1/mu0) E_n))",
    }

    # write JSON four levels up from verify/ -> .../LumenOS/Outputs/
    here = os.path.dirname(os.path.abspath(__file__))
    outputs = os.path.normpath(os.path.join(here, "..", "..", "..", "..", "Outputs"))
    os.makedirs(outputs, exist_ok=True)
    dest = os.path.join(outputs, "selection_rule_framing.json")
    with open(dest, "w") as f:
        json.dump(out, f, indent=2)

    # --- report ---
    print("=" * 64)
    print("A321 SELECTION-RULE FRAMING AUDIT (OI-287-1)")
    print("=" * 64)
    print(f"Taylor coefficients c1,c2,c3 = {c[0]:.6f}, {c[1]:.6f}, {c[2]:.6f}")
    print(f"  -> equal target (2,3,16): {out['taylor_match_target']}")
    print(f"Target normalized (c_n/c1): {target_norm.tolist()}  = (1, 1.5, 8)")
    print()
    print("Coefficient READOUT (corpus-faithful object):")
    print(f"  normalized = {(coeff_readout/coeff_readout[0]).tolist()}")
    print(f"  L2 to target = {out['coefficient_readout_L2_to_target']:.2e}  (TRIVIAL by definition)")
    print()
    print("Moments mu_n = int x^n rho dx (corpus's OTHER rho object, P20):")
    for n in range(4):
        print(f"  mu_{n} = {mu[n]:.6f}")
    print(f"  mu_0 = alpha^-1: {out['mu0_equals_alpha_inv']}")
    print(f"  ratios mu_n/mu_0 = {out['moment_ratios_over_mu0']}")
    print(f"  moment-ratio vec (n=1,2,3) normalized = {(mom_ratio_vec/mom_ratio_vec[0]).tolist()}")
    print(f"  -> do moments give (2,3,16)? {out['moments_give_target']}  "
          f"(L2={out['moments_L2_to_target']:.3f}, monotone-DECREASING)")
    print()
    print("Overlap-diagonal B[n] = <psi_n|c_k pi^k x^k|psi_n>/lambda_n (A296-A320 framing):")
    print(f"  eigenvalues = {[round(l,2) for l in lam]}")
    print(f"  normalized = {[round(b,3) for b in Bn]}")
    print(f"  mu entry = {Bn[1]:.3f} (target 1.5), L2 to target = {out['overlap_diagonal_L2_to_target']:.3f}")
    print()
    print("-" * 64)
    print(f"VERDICT: {out['verdict']}")
    print("-" * 64)
    print("The overlap-diagonal is A296's 'bare Taylor probe', EXCLUDED in")
    print("A296 (Third exclusion) and A297 (position moments), then revived")
    print("in A299 as a 'lead' and carried ~18 probes. The corpus mandates")
    print("NO overlap-diagonal: A287/A296 state the rule as a lepton MASS")
    print("correction law whose c_n are rho's Taylor coefficients, read")
    print("DIRECTLY. That readout gives (1,1.5,8) trivially. The magnitude")
    print("'gap' (mu=1.5 unreached across 7 arenas) is an artifact of")
    print("matching a functional the corpus excluded.")
    print()
    print(f"Wrote {dest}")
    return out


if __name__ == "__main__":
    main()
