verify-balanced-pairing.py

#!/usr/bin/env python3
"""Exact algebraic checks for the balanced-pairing normal-index bridge.

Standard library only. This checks identities, not topological instances.
The GHHR integral homology theorem and published boundary inversion identity
are explicit mathematical dependencies, not results proved by this program.
No manifolds, slopes, triangulations, or normal supports are enumerated.
"""
from collections import defaultdict
from fractions import Fraction as Q
from pathlib import Path
import argparse
import json


class Poly:
    """Rational commutative polynomials; adapted from the prior packet checker."""
    def __init__(self, terms=None):
        self.terms = {tuple(m): Q(c) for m, c in (terms or {}).items() if c}

    @staticmethod
    def of(x):
        return x if isinstance(x, Poly) else Poly({(): Q(x)})

    @staticmethod
    def var(name):
        return Poly({(name,): 1})

    def __add__(self, other):
        out = defaultdict(Q, self.terms)
        for m, c in Poly.of(other).terms.items():
            out[m] += c
        return Poly(out)

    __radd__ = __add__

    def __neg__(self):
        return Poly({m: -c for m, c in self.terms.items()})

    def __sub__(self, other):
        return self + (-Poly.of(other))

    def __rsub__(self, other):
        return Poly.of(other) - self

    def __mul__(self, other):
        out = defaultdict(Q)
        for m, c in self.terms.items():
            for n, d in Poly.of(other).terms.items():
                out[tuple(sorted(m+n))] += c*d
        return Poly(out)

    __rmul__ = __mul__

    def __truediv__(self, other):
        return self * (Q(1)/Q(other))

    def __eq__(self, other):
        return self.terms == Poly.of(other).terms


def P(raw):
    return tuple(v for j in range(0, len(raw), 3)
                 for v in (raw[j]-raw[j+1], raw[j+2]-raw[j+1]))


def add(a, b):
    return tuple(x+y for x, y in zip(a, b))


def scale(c, a):
    return tuple(c*x for x in a)


def omega(a, b):
    return sum(a[j+1]*b[j]-a[j]*b[j+1] for j in range(0, len(a), 2))


def reduce_raw(raw):
    maxima = tuple(max(raw[j:j+3]) for j in range(0, len(raw), 3))
    X = tuple(maxima[j//3]-raw[j] for j in range(len(raw)))
    return maxima, X


def parity_solve(rows, target):
    """GF(2) elimination: solve sum c_i rows_i = target, without searching."""
    width = len(target)
    if any(len(r) != width for r in rows):
        raise ValueError("row dimensions differ")
    basis = {}
    for idx, row in enumerate(rows):
        vec = sum((int(x) & 1) << j for j, x in enumerate(row))
        witness = 1 << idx
        while vec:
            p = vec.bit_length()-1
            if p in basis:
                v, w = basis[p]
                vec ^= v
                witness ^= w
            else:
                basis[p] = (vec, witness)
                break
    vec = sum((int(x) & 1) << j for j, x in enumerate(target))
    witness = 0
    while vec:
        p = vec.bit_length()-1
        if p not in basis:
            return {"supported": False, "rank_mod_2": len(basis)}
        v, w = basis[p]
        vec ^= v
        witness ^= w
    coefficients = [(witness >> i) & 1 for i in range(len(rows))]
    assert all((sum(c*row[j] for c, row in zip(coefficients, rows))-target[j]) % 2 == 0
               for j in range(width))
    return {"supported": True, "coefficients_mod_2": coefficients,
            "rank_mod_2": len(basis)}


def inv_poch(n, degree):
    out = [1]+[0]*degree
    for step in range(1, n+1):
        for d in range(step, degree+1):
            out[d] += out[d-step]
    return out


def index_m_zero(e, a_degree):
    """Certified coefficients for I_Delta(0,e) through A^a_degree.

    Only nonnegative n(n+1)/2 base exponents occur when m=0. The loop
    ends exactly when the next increasing base degree exceeds the cutoff.
    """
    out = defaultdict(int)
    n = max(0, -e)
    while 2*n*(n+1) <= a_degree:
        base = 2*n*(n+1)
        cutoff = (a_degree-base)//4
        a, b = inv_poch(n, cutoff), inv_poch(n+e, cutoff)
        for x in range(cutoff+1):
            for y in range(cutoff-x+1):
                out[base+4*(x+y)] += (-1)**n*a[x]*b[y]
        n += 1
    return {d: c for d, c in sorted(out.items()) if c}


def universal_checks():
    checks = {}

    def check(name, condition):
        if not condition:
            raise AssertionError(name)
        checks[name] = "pass"

    E, n, S, rho, r, T, delta = (Poly.var(x) for x in
        ("edge_sum", "meridian", "prime_sum", "longitude_prime_sum",
         "longitude_power", "maxima_sum", "double_arcs"))
    char_degree = 2*E+2*n-S
    complete_degree = char_degree-2*n
    chi_X = E-T
    M = (S+r*rho)/2-T
    raw_A_exponent = r*rho-complete_degree+r*n-2*M
    raw_i_exponent = r*rho-(2*E+2*n-S)-2*M
    normal_A_exponent = r*n-2*chi_X
    normal_i_exponent = 2*n-2*chi_X
    check("complete character has valuation 2E-S", complete_degree == 2*E-S)
    check("complete character phase is NOT just B^(2E-S)",
          (2*E+2*n-S)-(2*E-S) == 2*n)
    check("native all-orders A exponent matches normal weight",
          raw_A_exponent == normal_A_exponent)
    check("native i phase matches modulo four for integer meridian",
          raw_i_exponent-normal_i_exponent == -4*n)
    check("all-r normal summand degree", raw_A_exponent+2*delta
          == 2*(delta-chi_X)+r*n)
    check("dropping the Weyl exponent introduces an error",
          raw_A_exponent-r*n-normal_A_exponent == -r*n)
    check("dropping the central longitude scalar introduces an error",
          raw_A_exponent-r*rho-normal_A_exponent == -r*rho)

    # Symbolic surjectivity reconstruction from the integer E+T witness.
    x, y, z, d, l1, l2, l3 = (Poly.var(s) for s in
                             ("x", "y", "z", "tet_shift", "l1", "l2", "l3"))
    X = (x, y, z)
    Rl = (l1, l2, l3)
    raw_R = add(add(scale(-2, X), (d, d, d)), scale(-1, scale(r, Rl)))
    check("integer normal witness reconstructs exact lambda+rb",
          add(P(raw_R), scale(r, P(Rl))) == scale(-2, P(X)))
    check("surjectivity reconstruction uses only an integer edge witness",
          add(raw_R, scale(r, Rl)) == add(scale(-2, X), (d, d, d)))

    # One local coordinate fixture. It is not asserted to be a triangulation.
    X0 = (2, 0, 3, 0, 1, 2)
    lng = (0, 1, -1, 1, 0, 2)
    shifts = (1, 1, 1, 3, 3, 3)
    R0 = add(add(scale(-2, X0), shifts), scale(-1, lng))
    U0 = scale(Q(1, 2), add(R0, lng))
    maxima, Xout = reduce_raw(U0)
    check("half-integral raw representatives are kept exactly",
          maxima == (Q(1, 2), Q(3, 2)) and Xout == X0)
    check("root parity coset, not even lambda, is required",
          all(v % 2 == 0 for v in add(P(R0), P(lng)))
          and any(v % 2 for v in P(R0)))
    _, Xshifted = reduce_raw(add(U0, (Q(3, 2),)*3+(Q(-1, 2),)*3))
    check("normal class is invariant under fractional common shifts", Xshifted == Xout)
    check("two balanced inputs have identical support modulo two",
          all((a-b) % 2 == 0 for a, b in zip(add(P(R0), P(lng)),
                                           add(P(R0), scale(-1, P(lng))))))

    # Formal reindexing is checked WITHOUT assuming equality of local seeds.
    check("outer A exponent after r inversion and n reindexing",
          (-r)*(-n) == r*n)
    check("outer parity phase is unchanged modulo four",
          2*(-n)-2*n == -4*n)
    check("boundary vector after reindexing is its opposite",
          (r, n) == tuple(-s for s in (-r, -n)))
    tau = Poly.var("other_cusp_prime_sum")
    M_other = (S+tau+r*rho)/2-T
    check("other-cusp native scalar matches the same normal weight",
          tau+r*rho-complete_degree+r*n-2*M_other == normal_A_exponent)
    check("other-cusp i phase is exact modulo four",
          tau+r*rho-(2*E+2*n-S)-2*M_other-normal_i_exponent == -4*n)
    check("other-cusp reindexing and simultaneous inversion leave a partial reversal",
          tuple(-x for x in (-1, 1, n)) == (1, -1, -n))

    # Check the left-transpose character identity in bilinear symbols.
    vl, vu, lu = (Poly.var(s) for s in ("omega_v_lambda", "omega_v_u", "omega_lambda_u"))
    original_total = (vl+vu+lu)/2
    shifted_total = vl+vu/2+(-vl+lu)/2
    check("left seed translate and shifted average have identical Weyl phase",
          original_total == shifted_total)
    vw, wl, vl = (Poly.var(s) for s in ("omega_v_w", "omega_w_lambda", "omega_v_lambda"))
    check("lattice translation scalar is A^(-omega(v,w)/2)",
          (vl+wl)/2 == -vw/2+(vw+vl+wl)/2)
    check("same complete character requires exact orthogonality",
          (2*E+2*n-S)+(-2*n) == 2*E-S)

    # This fixture checks elimination only; no topological realization is claimed.
    parity_rows = [(1, 1, 0, 0), (0, 0, 2, 0)]
    supported = parity_solve(parity_rows, (1, 1, 0, 0))
    unsupported = parity_solve(parity_rows, (1, 0, 0, 0))
    check("GF2 support solver returns exact witness", supported["supported"])
    check("GF2 support solver rejects empty affine coset", not unsupported["supported"])

    positive = index_m_zero(1, 16)
    negative = index_m_zero(-1, 16)
    # In the present convention the term n=1 begins at A^4, not A^2:
    # I_Delta(0,-1) has leading -q = -A^4. The nonzero first negative
    # electric term is distinct from the leading one of I_Delta(0,1).
    check("seed is not invariant under root inversion",
          positive.get(0) == 1 and negative.get(0, 0) == 0
          and negative.get(4) == -1)
    check("seed inversion shortcut fails before any global summation",
          positive != negative)

    return {
        "scope": "universal identities and local algebra; no topological instances searched",
        "checks": checks,
        "exact_checks": len(checks),
        "normal_bridge": "F0(l^r)=sum_n (-1)^n A^(r*n) Z(-r*l-n*s)",
        "balanced_value": "0, using the proved normal bridge and published GHHR Remark 8.3",
        "same_complete_character_one_cusp": "all monomial translates kill the vacuum",
        "other_cusp_test": "(c_l/c_t)*sum_n (-1)^n A^n [Z(t-l-ns)-Z(-t-l-ns)], with nonvanishing unresolved",
        "external_mathematical_dependencies": [
            "GHHR Theorem 7.1: integral kernel and image of the normal-homology map",
            "GHHR Section 7.1: peripheral raw rows have mod-2 class zero",
            "GHHR Definition 8.1 and Remark 8.3: extended index and boundary inversion",
            "retained global all-input convergence and correct-sided translate theorem"
        ],
        "local_coordinate_fixture": {"X": X0, "R": R0, "longitude": lng,
            "U": list(map(str, U0)), "maxima": list(map(str, maxima)),
            "topological_realization_asserted": False},
        "seed_series_calibration": {
            "cutoff_A_degree": 16,
            "I_Delta_0_1": positive, "I_Delta_0_minus1": negative,
            "used_to_prove_global_boundary_inversion": False
        }
    }


def main():
    parser = argparse.ArgumentParser(description=__doc__)
    parser.add_argument("--parity", type=Path,
                        help="optional JSON with root rows and a target; performs GF(2) elimination only")
    args = parser.parse_args()
    result = universal_checks()
    if args.parity:
        data = json.loads(args.parity.read_text())
        result["supplied_parity_test"] = parity_solve(data["rows"], data["target"])
        result["supplied_parity_test"]["warning"] = "not a triangulation or character-consistency verifier"
    output = json.dumps(result, indent=2, sort_keys=True)+"\n"
    (Path(__file__).resolve().parent/"exact-results.json").write_text(output)
    print(output, end="")


if __name__ == "__main__":
    main()