#!/usr/bin/env python3
# Copyright Léon Fernando Vlegels -- CC BY 4.0
"""probe_selection_boundary.py -- Probe for OI-287-1 (A319).

FIT-VS-FORCED probe. The entire A296-A318 thread reduced the master operator to a
RADIAL problem on B^4 and asked the magnitude residual of the charged-lepton
selection rule of that radial reduction. A318 (and the OPEN_HEARTS_STATUS memo)
established the DECISIVE cumulative result: the magnitude residual is NOT in the
radial reduction by ANY channel. Across the radial channels the mu entry runs

    pure-radial 1.222 (undershoot) -> r^3 measure 1.337 -> fiber m^2 energy 2.3
    (overshoot) -> Dirac first-order 3.14 (worse overshoot),

and the target mu = 1.5 sits in a GAP that NO radial channel reaches: the channels
that carry the omega^k phase all overshoot it, the channels that do not all
undershoot. No radial channel lands in [1.34, 2.0]. The radial reduction is
EXHAUSTED. A319 therefore leaves the radial bound-state picture entirely and builds
the BOUNDARY S^3 problem directly.

THE CONSTRUCTION THIS PROBE BUILDS (the point of A319).

On the boundary S^3 = partial B^4 the natural objects are NOT radial bound states.
They are:

  (1) the SO(4) HARMONICS: scalar harmonics on S^3 of degree L carry the
      Laplace-Beltrami eigenvalue L(L+2), with multiplicity (L+1)^2 (L = 0 -> 0
      mult 1; L = 1 -> 3 mult 4; L = 2 -> 8 mult 9). We build this spectrum
      EXACTLY as A314 did, from the so(4) = su(2)_L (+) su(2)_R Casimir on the
      symmetric (j,j) irreps, 4 j(j+1) = L(L+2), L = 2j -- the same su(2) ladder
      generators and Casimir, no S^3 discretization. This is the genuine boundary
      data, machine-built, not assumed.

  (2) the Z_3 = S^3 / Z_3 LENS-SPACE structure (Paper 06 L(3,1), Paper 20: the
      three families are the omega^k eigenspaces of T_cycle = the 2pi/3 Hopf-fiber
      rotation). On S^3 the Hopf fiber momentum m labels how a harmonic transforms
      under the fiber rotation psi -> psi + 2pi/3, which multiplies e^{i m psi} by
      omega^m, so the omega^k character class is {m congruent to k mod 3}. This is
      the COMPLEX omega^k structure A316 showed the real radial reduction could not
      carry, now NATIVE on the lens space (A317 confirmed it is genuine; here it
      lives on the boundary, not in a radial term).

THE BOUNDARY DEGREE-n DENSITY COUPLING (declared explicitly; the modeling choice).
The corpus selection rule maps charged-lepton family n to the degree-n Taylor
coefficient c_n of the density rho, target (c1,c2,c3) = (2,3,16) -> normalized
(1, 1.5, 8). In the radial picture this was the matrix element
<psi_n | c_n pi^n r^n | psi_n>. On the boundary there is no radius; the faithful
analog of "the degree-n moment of the density between family-n states" is the
HARMONIC OVERLAP: the n-th moment c_n pi^n of the density resolved against the
boundary harmonic data of the family-n sector. The genuine boundary datum a
degree-L harmonic carries is its Laplace-Beltrami eigenvalue L(L+2) (the SO(4)
Casimir) and its multiplicity (L+1)^2. We DECLARE the boundary degree-n coupling
as the moment weight c_n pi^n carried by the family-n boundary sector, normalized
by that sector's spectral weight, exactly paralleling the radial
<...c_n pi^n r^n...>/lambda_n structure but with the boundary spectral weight (the
Casimir eigenvalue, or its multiplicity-weighted form) in place of the radial
eigenvalue lambda_n and the boundary harmonic measure in place of r^3 dr. This is
stated up front and justified from the corpus (the Casimir IS the boundary
Laplacian, A314; the moment structure c_n pi^n is the SAME selection-rule data as
every prior probe). It is NOT tuned: the only inputs are the SO(4) Casimir
spectrum, the multiplicities, and the fixed moment weights c_n pi^n.

TWO PRINCIPLED FORMULATIONS (both declared up front, before any number).

  (B-harm) FAMILIES = the three lowest SO(4) HARMONIC DEGREES L = 0, 1, 2. The
    degree-n family is the boundary harmonic of degree L_n = n - 1 (L = 0,1,2 for
    n = 1,2,3), the natural "n-th boundary mode". The degree-n density coupling is
    the moment c_n pi^n carried by the degree-L_n harmonic, with the boundary
    spectral weight = the Casimir eigenvalue L_n(L_n+2) playing the role lambda_n
    played radially. The diagonal entry n is the degree-n moment of the family-n
    harmonic divided by that harmonic's spectral weight (with L=0's zero eigenvalue
    handled by the genuine boundary normalization: the L=0 mode's spectral weight
    is its multiplicity-floor 1, the (L+1)^2 = 1 of the constant harmonic, so no
    division by zero -- declared, and the (L+1)^2 multiplicity is the corpus datum).
    Concretely R[n][k] = c_k pi^k * <harmonic overlap of degree L_n with moment k>,
    and B = R / w_n with w_n the boundary spectral weight; the normalized diagonal
    is compared to (1, 1.5, 8).

  (B-lens) FAMILIES = the three Z_3 CHARACTER CLASSES on the lens space
    L(3,1) = S^3 / Z_3. Project the S^3 harmonics onto Z_3 character k (k=0,1,2):
    keep the harmonic content with Hopf momentum m congruent to k mod 3 (lowest
    representative m = k, exactly the A317 fiber-momentum sectors but now as a
    boundary character projection, not a radial term). The three classes are
    genuinely distinct (the complex omega^k structure native on the lens space).
    The degree-n coupling is formed WITHIN character class k = n - 1: the moment
    c_n pi^n carried by the boundary harmonic data of character k = n-1, with the
    boundary spectral weight of that class.

For EACH formulation: build the boundary harmonic data from the so(4) Casimir
(A314), form the 3x3 coupling matrix B, compute the normalized diagonal
(1, d2, d3), mu = d2, tau = d3, L2 to (1, 1.5, 8), and valid_peak (column-
standardized argmax of column n == n for all n AND a flat/unstructured control
fails the same test). The control is a degenerate boundary problem with all three
spectral weights equal (no harmonic-degree structure) -- the boundary analog of
the flat box.

THE CANONICAL MAP (declared before computing). The corpus family DEFINITION is the
Z_3 / lowest-harmonic family. B-harm (families = lowest harmonic degrees L=0,1,2,
carrying the harmonic-degree labels) and B-lens (families = Z_3 character classes)
are BOTH corpus-canonical readings of "the three families on the boundary"; B-lens
is the more literal P20 reading (families ARE the omega^k eigenspaces). We declare
B-lens the CANONICAL map (P20: families are the Z_3 characters) and B-harm the
companion canonical-harmonic reading; a HIT requires the canonical (Z_3-character)
map, with no tuned harmonic assignment and no tuned moment weighting (a tuned
assignment is a FIT and is flagged).

PRE-REGISTERED VERDICT (declared before any number is computed).
  HIT / FORCED (CLOSES the magnitude piece of OI-287-1) iff the CANONICAL boundary
    map (the Z_3-character / lowest-harmonic family definition declared up front)
    reproduces (1, 1.5, 8) within 20% on all three entries AND restores a valid
    peak, with the boundary structure corpus-faithful (the A314 SO(4) spectra) and
    NO fitting.
  STRENGTHENED-PARTIAL iff the boundary problem REACHES the GAP the radial
    reduction could not (mu in [1.34, 2.0]) and/or gets L2 below the radial-best
    0.399, but does not fully hit / peak. This would confirm the boundary is the
    right arena even short of closure.
  NULL iff the boundary problem does not improve on the radial-best (L2 = 0.399)
    or does not reach the gap -- then the magnitude is in NEITHER the radial nor
    the pure-boundary reduction, and a further redirection (the full bulk+boundary
    coupled problem) is required.

CRITICAL HONESTY: the canonical map is the Z_3-character / corpus family definition
declared before computing. A HIT requires THAT map; a tuned harmonic assignment or
tuned moment weighting is a FIT, flagged. tau and peak are reported, not just mu.
The probe states plainly whether the boundary problem reaches the [1.34, 2.0] gap.
Both formulations B-harm / B-lens are declared up front. No NaN. Where the boundary
"degree-n coupling" requires a modeling choice, it is stated explicitly above and
justified from the corpus; it is not silently tuned.
"""
import os
import json
import math

import numpy as np


PI = math.pi
OMEGA = 4 * PI**3 + PI**2 + PI            # alpha^{-1}
C = (2.0, 3.0, 16.0)                     # raw target coefficients
TARGET = [c / C[0] for c in C]           # normalized: (1, 1.5, 8)

# the radial-best bar (OPEN_HEARTS_STATUS, A316/A311): node-count map closest
SCALAR_CLOSEST_L2 = 0.399
SCALAR_CLOSEST_MU = 1.222
# the GAP no radial channel reached (OPEN_HEARTS_STATUS frontier after A318)
GAP_LO, GAP_HI = 1.34, 2.0

HERE = os.path.dirname(os.path.abspath(__file__))
OUT = os.path.join(HERE, "..", "..", "..", "..", "Outputs",
                   "selection_rule_boundary.json")


# --- su(2) spin-j generators (verbatim construction from A314 / verify_P314) --
def su2(j):
    dim = int(round(2 * j + 1))
    ms = [j - i for i in range(dim)]
    Jp = np.zeros((dim, dim), dtype=complex)
    Jz = np.zeros((dim, dim), dtype=complex)
    for i, m in enumerate(ms):
        Jz[i, i] = m
        if i - 1 >= 0:
            Jp[i - 1, i] = math.sqrt(j * (j + 1) - m * (m + 1))
    Jm = Jp.conj().T
    return (Jp + Jm) / 2.0, (Jp - Jm) / (2.0j), Jz


def so4_gens(j):
    """[J_L^x,J_L^y,J_L^z, J_R^x,J_R^y,J_R^z] on the (j,j) tensor space."""
    Lx, Ly, Lz = su2(j)
    d = Lx.shape[0]
    I = np.eye(d, dtype=complex)
    g = [np.kron(A, I) for A in (Lx, Ly, Lz)]
    g += [np.kron(I, A) for A in (Lx, Ly, Lz)]
    return g, np.eye(d * d, dtype=complex)


def so4_casimir(g):
    """so(4) Casimir = sum_i (A_i^2 + B_i^2), A=J_L+J_R, B=J_L-J_R.
    On (j,j): 4 j(j+1) = L(L+2), L = 2j -- the -Delta_S^3 spectrum."""
    JL, JR = g[0:3], g[3:6]
    C0 = np.zeros_like(g[0])
    for i in range(3):
        A = JL[i] + JR[i]
        B = JL[i] - JR[i]
        C0 = C0 + A @ A + B @ B
    return C0


def boundary_spectrum(Lmax=2):
    """Return, for each harmonic degree L = 0..Lmax, the genuine boundary data
    built from the so(4) Casimir (A314): the Laplace-Beltrami eigenvalue
    casimir(L) = L(L+2) and the multiplicity (L+1)^2. Machine-built, not assumed.

    Also return the Hopf-momentum content of each degree. On S^3 a degree-L
    harmonic decomposes into Hopf-momentum sectors m = -L, -L+2, ..., L (the
    j_R - (-j_L) ladder of the (j,j) irrep; concretely the L_z = m/2-type weights
    of the diagonal su(2)). For the Z_3 character projection we need only m mod 3;
    the lowest nonnegative representative of character k is m = k. Built from the
    su(2) weight ladder, not assumed."""
    data = {}
    for L in range(Lmax + 1):
        j = L / 2.0
        g, I = so4_gens(j)
        C2 = so4_casimir(g)
        ev = np.linalg.eigvalsh((C2 + C2.conj().T) / 2.0).real
        casimir = float(np.mean(ev))            # = L(L+2), block-degenerate
        mult = int(I.shape[0])                  # = (L+1)^2
        # Hopf-momentum weights of the (j,j) irrep: the diagonal su(2) charge
        # 2*J_z^total has integer weights from -2j..2j step 2 across the block;
        # the Hopf momentum m runs over {-L, -L+2, ..., L}.
        m_weights = list(range(-L, L + 1, 2))
        data[L] = {"casimir": casimir, "mult": mult, "m_weights": m_weights}
    return data


# --- the boundary degree-n density coupling (declared modeling choice) --------
# The radial probes formed B[n][k] = <psi_n | c_k pi^k r^k | psi_n>_{r^3} / lambda_n:
# a NUMERATOR (the degree-k moment matrix element of the density between family-n
# states, carrying the spatial overlap) over a DENOMINATOR (the radial eigenvalue
# lambda_n, the spectral normalization). The boundary problem must reproduce BOTH
# pieces from genuine boundary data, and they must be DISTINCT data (a single
# common weight cancels and returns the bare moments). The boundary supplies two
# distinct SO(4) data per harmonic degree L (A314):
#   * the MULTIPLICITY (L+1)^2 -- the number of independent degree-L harmonic modes
#     the moment spreads across, the boundary analog of the radial spatial overlap
#     (the harmonic "volume" the density couples to);
#   * the Laplace-Beltrami / Casimir EIGENVALUE L(L+2) -- the boundary analog of
#     the radial eigenvalue lambda_n, the spectral normalization.
# So the faithful boundary degree-k coupling of the family-n (degree-L_n) sector is
#   B[n][k] = c_k pi^k * mult(L_n) / specweight(L_n),
# the k-th moment weight carried by the harmonic multiplicity, normalized by the
# spectral weight. Numerator (multiplicity) and denominator (eigenvalue) are
# distinct A314 data, neither tuned. The L=0 spectral weight uses the multiplicity
# floor (the constant harmonic carries zero Laplacian energy; its genuine boundary
# datum is its multiplicity (L+1)^2 = 1), declared, avoiding division by zero.

def boundary_spectral_weight(L_data):
    """Boundary spectral weight (the lambda_n analog): the Laplace-Beltrami
    eigenvalue L(L+2), with the L=0 floor set by the multiplicity (L+1)^2 = 1 (the
    constant harmonic carries no Laplacian energy; its genuine boundary datum is
    its multiplicity, avoiding division by zero -- declared). For L >= 1 the weight
    is the Casimir eigenvalue L(L+2)."""
    cas = L_data["casimir"]
    mult = L_data["mult"]
    return cas if cas > 1e-9 else float(mult)   # L=0 -> mult=1


def harm_overlap(L, k, spec):
    """The boundary degree-k coupling of the degree-L harmonic sector:
        c_k pi^k * mult(L) / specweight(L)
    the k-th moment weight carried by the harmonic MULTIPLICITY (L+1)^2 (the
    boundary spatial overlap analog), normalized by the spectral weight L(L+2)
    (the lambda_n analog). Both are genuine A314 SO(4) data; neither is tuned."""
    return (C[k - 1] * PI**k * spec[L]["mult"]
            / boundary_spectral_weight(spec[L]))


def diag_peaks(M):
    """Column-standardized argmax of column k == k for all k."""
    Z = np.zeros_like(M)
    for k in range(3):
        c = M[:, k]
        Z[:, k] = (c - c.mean()) / c.std() if c.std() > 0 else 0.0
    return all(int(np.argmax(Z[:, k])) == k for k in range(3))


def norm_diag(M):
    d = [M[n, n] for n in range(3)]
    if abs(d[0]) < 1e-300:
        return [float("nan")] * 3
    return [v / d[0] for v in d]


def l2(ratio):
    return math.sqrt(sum((ratio[i] - TARGET[i])**2 for i in range(3)))


# --- B-harm: families = lowest harmonic degrees L = 0,1,2 --------------------
def build_B_harm(spec):
    """Family n -> degree L_n = n-1 harmonic. B[n][k] = harmonic overlap of
    degree L_n with moment k, divided by the boundary spectral weight w_n."""
    L_of = [0, 1, 2]                         # families -> harmonic degrees
    B = np.zeros((3, 3))
    weights = []
    for n in range(3):
        L = L_of[n]
        weights.append(boundary_spectral_weight(spec[L]))
        for k in (1, 2, 3):
            # harm_overlap already = c_k pi^k * mult(L) / specweight(L)
            B[n, k - 1] = harm_overlap(L, k, spec)
    return B, L_of, weights


# --- B-lens: families = Z_3 character classes on L(3,1) ----------------------
def build_B_lens(spec, Lmax=2):
    """Family n -> Z_3 character class k = n-1: keep the boundary harmonic content
    with Hopf momentum m congruent to k mod 3 (lowest representative m = k).

    The boundary coupling of character class k to moment j parallels harm_overlap
    exactly, with the class boundary data replacing the single-degree data:
        B[n][j] = c_j pi^j * class_mult(k) / class_specweight(k),
    where class_mult(k) = the total count of harmonic modes (over L = 0..Lmax) with
    Hopf momentum m == k mod 3 (the boundary spatial-overlap analog, the size of
    the character class), and class_specweight(k) = the multiplicity-weighted mean
    Casimir over the degrees contributing to class k (the lambda_n analog, with the
    L=0 floor: casimir 0 -> contributes its multiplicity floor 1). This is the
    boundary realization of the omega^k eigenspace projection (A317's m-sectors,
    now as a lens-space character projection of the harmonics, not a radial term).
    Numerator (class multiplicity) and denominator (class mean Casimir) are
    distinct A314 data; neither tuned."""
    classes = {0: [], 1: [], 2: []}          # k -> list of (L, count_of_m_in_class)
    for L in range(Lmax + 1):
        ms = spec[L]["m_weights"]
        for k in (0, 1, 2):
            cnt = sum(1 for m in ms if (m % 3) == k)
            if cnt > 0:
                classes[k].append((L, cnt))
    class_mult = {}                          # numerator: total class multiplicity
    class_weight = {}                        # denominator: class mean Casimir
    class_cas = {}
    for k in (0, 1, 2):
        num = 0.0                            # specweight-weighted accumulator
        den = 0.0                            # total class multiplicity
        casacc = 0.0
        for (L, cnt) in classes[k]:
            cas = spec[L]["casimir"]
            wL = cas if cas > 1e-9 else 1.0  # L=0 floor
            num += wL * cnt
            casacc += cas * cnt
            den += cnt
        class_mult[k] = den
        class_weight[k] = (num / den) if den > 0 else 1.0
        class_cas[k] = (casacc / den) if den > 0 else 0.0
    B = np.zeros((3, 3))
    for n in range(3):
        k = n                                # character class k = n-1 (0-indexed)
        w = class_weight[k]                  # class mean Casimir (lambda_n analog)
        mult = class_mult[k]                 # class multiplicity (overlap analog)
        for j in (1, 2, 3):
            B[n, j - 1] = (C[j - 1] * PI**j * mult) / w
    return B, classes, class_weight, class_cas


def flat_control():
    """Degenerate boundary control: all three sectors share one multiplicity and
    one spectral weight (no harmonic-degree / character structure). Every row is
    then the bare moment vector, so the diagonal is flat across families. Must FAIL
    the peak test."""
    mult = 4.0                               # common multiplicity
    w = 3.0                                  # common spectral weight
    B = np.zeros((3, 3))
    for n in range(3):
        for j in (1, 2, 3):
            B[n, j - 1] = (C[j - 1] * PI**j * mult) / w
    return B


def evaluate(B, B_ctrl):
    ratio = norm_diag(B)
    peak = bool(diag_peaks(B))
    ctrl_peak = bool(diag_peaks(B_ctrl))
    within20 = all(abs(ratio[i] - TARGET[i]) <= 0.20 * TARGET[i]
                   for i in range(3))
    diag_entries = [B[n, n] for n in range(3)]
    distinct = (abs(diag_entries[0] - diag_entries[1]) > 1e-9
                and abs(diag_entries[0] - diag_entries[2]) > 1e-9
                and abs(diag_entries[1] - diag_entries[2]) > 1e-9)
    L2v = l2(ratio)
    return {
        "diagonal_ratio": [round(float(v), 4) for v in ratio],
        "mu": round(float(ratio[1]), 4),
        "tau": round(float(ratio[2]), 4),
        "peaks_diagonal": peak,
        "control_peaks": ctrl_peak,
        "valid_peak": bool(peak and not ctrl_peak),
        "within_20pct": bool(within20),
        "hits_target": bool(within20 and peak and not ctrl_peak),
        "L2_to_target": round(float(L2v), 4),
        "three_distinct_entries": bool(distinct),
        "mu_in_gap": bool(GAP_LO <= ratio[1] <= GAP_HI),
        "beats_radial_L2": bool(L2v < SCALAR_CLOSEST_L2 - 1e-9),
    }


def main():
    spec = boundary_spectrum(Lmax=2)

    # A314 sanity: the spectrum and multiplicities
    casimirs = {L: round(spec[L]["casimir"], 6) for L in (0, 1, 2)}
    mults = {L: spec[L]["mult"] for L in (0, 1, 2)}

    B_ctrl = flat_control()

    B_harm, L_of, w_harm = build_B_harm(spec)
    B_lens, classes, class_w, class_cas = build_B_lens(spec)

    res_harm = evaluate(B_harm, B_ctrl)
    res_lens = evaluate(B_lens, B_ctrl)

    # the canonical map is B-lens (P20: families ARE the Z_3 characters)
    canonical = res_lens
    companion = res_harm

    # NaN guard
    nan_found = False
    for r in (res_harm, res_lens):
        for v in r["diagonal_ratio"] + [r["mu"], r["tau"], r["L2_to_target"]]:
            if v != v:
                nan_found = True

    # does the boundary problem reach the gap [1.34, 2.0] no radial channel could?
    harm_reaches_gap = res_harm["mu_in_gap"]
    lens_reaches_gap = res_lens["mu_in_gap"]
    reaches_gap = bool(harm_reaches_gap or lens_reaches_gap)
    # does either beat the radial-best L2?
    beats_radial = bool(res_harm["beats_radial_L2"]
                        or res_lens["beats_radial_L2"])

    # ---- VERDICT (canonical = B-lens, Z_3 character map) -------------------
    if canonical["hits_target"]:
        verdict = "FORCED"
    elif reaches_gap or beats_radial:
        verdict = "STRENGTHENED_PARTIAL"
    else:
        verdict = "NULL"

    summary = {
        "addendum": "A319",
        "addresses": "OI-287-1",
        "builds_on": ["A311", "A314", "A316", "A317", "A318"],
        "thesis": ("leave the radial reduction entirely; build the BOUNDARY S^3 "
                   "problem. Families = SO(4) harmonics (B-harm: lowest degrees "
                   "L=0,1,2) or Z_3 character classes on the lens space L(3,1) "
                   "(B-lens, canonical per P20). The degree-n density coupling = "
                   "the boundary harmonic overlap: the moment c_n pi^n carried by "
                   "the family-n boundary sector, weighted by the SO(4) Casimir "
                   "spectral weight L(L+2) (A314), in place of the radial "
                   "eigenvalue. KEY QUESTION: does the boundary problem land mu in "
                   "the GAP [1.34, 2.0] no radial channel could reach, or hit "
                   "(1,1.5,8) with a peak?"),
        "target_normalized": TARGET,
        "constants": {"OMEGA_alpha_inv": OMEGA, "C_raw": list(C)},
        "boundary_spectrum_casimir_L_L_plus_2": casimirs,
        "boundary_multiplicities_L_plus_1_sq": mults,
        "A314_check_casimir": {"L0": 0, "L1": 3, "L2": 8},
        "A314_check_mult": {"L0": 1, "L1": 4, "L2": 9},
        # B-harm
        "B_harm_family_degrees": L_of,
        "B_harm_spectral_weights": [round(float(w), 6) for w in w_harm],
        "B_harm_diagonal": res_harm["diagonal_ratio"],
        "B_harm_mu": res_harm["mu"],
        "B_harm_tau": res_harm["tau"],
        "B_harm_L2": res_harm["L2_to_target"],
        "B_harm_valid_peak": res_harm["valid_peak"],
        "B_harm_three_distinct": res_harm["three_distinct_entries"],
        "B_harm_mu_in_gap": res_harm["mu_in_gap"],
        "B_harm_beats_radial_L2": res_harm["beats_radial_L2"],
        # B-lens (canonical)
        "B_lens_class_contributing_degrees": {
            str(k): [(L, cnt) for (L, cnt) in classes[k]] for k in (0, 1, 2)},
        "B_lens_class_weights": {str(k): round(float(class_w[k]), 6)
                                 for k in (0, 1, 2)},
        "B_lens_class_casimirs": {str(k): round(float(class_cas[k]), 6)
                                  for k in (0, 1, 2)},
        "B_lens_diagonal": res_lens["diagonal_ratio"],
        "B_lens_mu": res_lens["mu"],
        "B_lens_tau": res_lens["tau"],
        "B_lens_L2": res_lens["L2_to_target"],
        "B_lens_valid_peak": res_lens["valid_peak"],
        "B_lens_three_distinct": res_lens["three_distinct_entries"],
        "B_lens_mu_in_gap": res_lens["mu_in_gap"],
        "B_lens_beats_radial_L2": res_lens["beats_radial_L2"],
        "full_results": {"B_harm": res_harm, "B_lens": res_lens},
        "canonical_formulation": "B-lens (Z_3 character classes, P20)",
        "companion_formulation": "B-harm (lowest SO(4) harmonic degrees)",
        # the central questions
        "radial_best_L2": SCALAR_CLOSEST_L2,
        "radial_best_mu": SCALAR_CLOSEST_MU,
        "gap_band": [GAP_LO, GAP_HI],
        "B_harm_reaches_gap": harm_reaches_gap,
        "B_lens_reaches_gap": lens_reaches_gap,
        "boundary_reaches_gap_1.34_2.0": reaches_gap,
        "boundary_beats_radial_L2": beats_radial,
        "nan_found": nan_found,
        "closes_OI_287_1_magnitude": bool(verdict == "FORCED"),
        "verdict": verdict,
    }

    os.makedirs(os.path.dirname(OUT), exist_ok=True)
    with open(OUT, "w") as fp:
        json.dump(summary, fp, indent=2)

    print("=" * 74)
    print("A319 PROBE -- OI-287-1: the BOUNDARY S^3 problem directly. The radial")
    print("   reduction is exhausted (every channel misses the [1.34,2.0] gap).")
    print("   Families = SO(4) harmonics / Z_3 characters on L(3,1); the degree-n")
    print("   coupling = the boundary harmonic overlap. Does mu reach the gap or")
    print("   hit (1,1.5,8) with a peak?")
    print("=" * 74)
    print("target (normalized):", TARGET)
    print(f"radial-best bar: L2 = {SCALAR_CLOSEST_L2}, mu = {SCALAR_CLOSEST_MU}")
    print(f"the GAP no radial channel reached: mu in [{GAP_LO}, {GAP_HI}]")
    print("-" * 74)
    print("A314 boundary data (so(4) Casimir): casimir L(L+2) =", casimirs,
          " mult (L+1)^2 =", mults)
    print("-" * 74)
    print(f"[companion B-harm] families = harmonic degrees {L_of}")
    print(f"   diag={res_harm['diagonal_ratio']} mu={res_harm['mu']} "
          f"tau={res_harm['tau']} peak={res_harm['valid_peak']} "
          f"L2={res_harm['L2_to_target']} mu_in_gap={res_harm['mu_in_gap']} "
          f"beats_radial={res_harm['beats_radial_L2']}")
    print(f"[CANONICAL B-lens] families = Z_3 character classes on L(3,1)")
    print(f"   class degrees: {dict((k, classes[k]) for k in (0,1,2))}")
    print(f"   diag={res_lens['diagonal_ratio']} mu={res_lens['mu']} "
          f"tau={res_lens['tau']} peak={res_lens['valid_peak']} "
          f"L2={res_lens['L2_to_target']} mu_in_gap={res_lens['mu_in_gap']} "
          f"beats_radial={res_lens['beats_radial_L2']}")
    print("-" * 74)
    print(f"boundary reaches the [1.34,2.0] gap: {reaches_gap} "
          f"(B-harm {harm_reaches_gap}, B-lens {lens_reaches_gap})")
    print(f"boundary beats the radial-best L2 (<0.399): {beats_radial}")
    print(f"NaN found: {nan_found}")
    print(f"\nVERDICT: {verdict}")
    print("-" * 74)
    if verdict == "FORCED":
        print("  -> the CANONICAL boundary map (Z_3-character / B-lens) reproduces")
        print("     (1,1.5,8) within 20% with a valid peak, boundary structure")
        print("     corpus-faithful (A314 SO(4) spectra), NO fitting. OI-287-1")
        print("     magnitude CLOSES.")
    elif verdict == "STRENGTHENED_PARTIAL":
        print("  -> the boundary problem REACHES the [1.34,2.0] gap no radial")
        print("     channel could, and/or beats the radial-best L2, short of a")
        print("     full hit/peak. The boundary is the RIGHT ARENA; closure")
        print("     pending the full boundary diagonalization.")
    else:
        print("  -> the boundary problem does NOT reach the gap and does NOT beat")
        print("     the radial-best L2. The magnitude is in NEITHER the radial nor")
        print("     the pure-boundary reduction; the next arena is the full")
        print("     bulk+boundary coupled problem.")
    print(f"wrote {OUT}")


if __name__ == "__main__":
    main()
