verify-peripheral-filling-base-change.py

#!/usr/bin/env python3
"""Exact checks for the empty-sector peripheral filling ideal.

The script uses only integer arithmetic.  It verifies:

1. the Frohman--Gelca product-to-sum law in the rank-two quantum torus;
2. covariance of the two solid-torus filling generators under SL(2,Z);
3. the exact classical factorization showing that the second generator is
   radical-redundant at t=1, although it is essential quantum mechanically.

These identities do not prove ordinary-skein excision or a filled quantum
trace.  The companion obstruction note explains why empty-sector peripheral
relations are insufficient.  No numerical approximation or manifold census
is used.
"""

from __future__ import annotations

from dataclasses import dataclass
from typing import Dict, Iterable, Tuple


Vector = Tuple[int, int]
QKey = Tuple[int, int, int]  # t exponent, first torus exponent, second torus exponent
LKey = Tuple[int, int]       # commuting eigenvalue exponents


def det(u: Vector, v: Vector) -> int:
    return u[0] * v[1] - u[1] * v[0]


@dataclass(frozen=True)
class QuantumElement:
    """Finite Z-linear combinations of t^k E_(a,b).

    Multiplication is E_u E_v = t^det(u,v) E_(u+v).
    """

    terms: Dict[QKey, int]

    def __post_init__(self) -> None:
        clean = {key: value for key, value in self.terms.items() if value}
        object.__setattr__(self, "terms", clean)

    def __add__(self, other: "QuantumElement") -> "QuantumElement":
        result = dict(self.terms)
        for key, value in other.terms.items():
            result[key] = result.get(key, 0) + value
        return QuantumElement(result)

    def __sub__(self, other: "QuantumElement") -> "QuantumElement":
        return self + other.scale(-1)

    def scale(self, scalar: int) -> "QuantumElement":
        return QuantumElement({key: scalar * value for key, value in self.terms.items()})

    def __mul__(self, other: "QuantumElement") -> "QuantumElement":
        result: Dict[QKey, int] = {}
        for (kt, a, b), c in self.terms.items():
            for (lt, x, y), d in other.terms.items():
                exponent = kt + lt + det((a, b), (x, y))
                key = (exponent, a + x, b + y)
                result[key] = result.get(key, 0) + c * d
        return QuantumElement(result)

    def apply_sl2(self, matrix: Tuple[Tuple[int, int], Tuple[int, int]]) -> "QuantumElement":
        (a, b), (c, d) = matrix
        if a * d - b * c != 1:
            raise ValueError("matrix is not in SL(2,Z)")
        return QuantumElement(
            {
                (kt, a * x + b * y, c * x + d * y): coefficient
                for (kt, x, y), coefficient in self.terms.items()
            }
        )


def e(u: Vector, t_power: int = 0) -> QuantumElement:
    return QuantumElement({(t_power, u[0], u[1]): 1})


def scalar_t(t_power: int, coefficient: int = 1) -> QuantumElement:
    return QuantumElement({(t_power, 0, 0): coefficient})


def k(u: Vector, t_power: int = 0) -> QuantumElement:
    return e(u, t_power) + e((-u[0], -u[1]), t_power)


def add_vectors(u: Vector, v: Vector) -> Vector:
    return (u[0] + v[0], u[1] + v[1])


def subtract_vectors(u: Vector, v: Vector) -> Vector:
    return (u[0] - v[0], u[1] - v[1])


def product_to_sum_rhs(u: Vector, v: Vector) -> QuantumElement:
    intersection = det(u, v)
    return k(add_vectors(u, v), intersection) + k(subtract_vectors(u, v), -intersection)


def filling_generators(rho: Vector, slope: Vector) -> Tuple[QuantumElement, QuantumElement]:
    """Frohman--Gelca generators in their t-convention.

    The ordered basis satisfies det(rho, slope)=1.  The first relation is the
    meridian disk relation; the second records the framing change.
    """

    if det(rho, slope) != 1:
        raise ValueError("expected det(rho, slope)=1")
    first = k(slope) + scalar_t(2) + scalar_t(-2)
    second = k(add_vectors(rho, slope)) + k(rho, -3)
    return first, second


@dataclass(frozen=True)
class LaurentElement:
    terms: Dict[LKey, int]

    def __post_init__(self) -> None:
        clean = {key: value for key, value in self.terms.items() if value}
        object.__setattr__(self, "terms", clean)

    def __add__(self, other: "LaurentElement") -> "LaurentElement":
        result = dict(self.terms)
        for key, value in other.terms.items():
            result[key] = result.get(key, 0) + value
        return LaurentElement(result)

    def __mul__(self, other: "LaurentElement") -> "LaurentElement":
        result: Dict[LKey, int] = {}
        for (a, b), c in self.terms.items():
            for (x, y), d in other.terms.items():
                key = (a + x, b + y)
                result[key] = result.get(key, 0) + c * d
        return LaurentElement(result)


def monomial(a: int, b: int, coefficient: int = 1) -> LaurentElement:
    return LaurentElement({(a, b): coefficient})


def classical_k(u: Vector) -> LaurentElement:
    return monomial(*u) + monomial(-u[0], -u[1])


def exact_product_to_sum_checks() -> int:
    checked = 0
    vectors: Iterable[Vector] = (
        (-3, 2), (-2, -1), (-1, 3), (0, 1), (1, 0), (1, 2), (2, -3), (3, 1)
    )
    for u in vectors:
        for v in vectors:
            assert k(u) * k(v) == product_to_sum_rhs(u, v)
            checked += 1
    return checked


def exact_sl2_checks() -> int:
    matrices = (
        ((1, 0), (0, 1)),
        ((0, -1), (1, 0)),
        ((1, 1), (0, 1)),
        ((2, 1), (1, 1)),
        ((-1, 0), (0, -1)),
    )
    rho, slope = (1, 0), (0, 1)
    standard = filling_generators(rho, slope)
    checked = 0
    for matrix in matrices:
        transformed_rho = (
            matrix[0][0] * rho[0] + matrix[0][1] * rho[1],
            matrix[1][0] * rho[0] + matrix[1][1] * rho[1],
        )
        transformed_slope = (
            matrix[0][0] * slope[0] + matrix[0][1] * slope[1],
            matrix[1][0] * slope[0] + matrix[1][1] * slope[1],
        )
        transported = tuple(generator.apply_sl2(matrix) for generator in standard)
        direct = filling_generators(transformed_rho, transformed_slope)
        assert transported == direct
        checked += 1
    return checked


def exact_classical_factorization_checks() -> Tuple[LaurentElement, LaurentElement]:
    rho, slope = (1, 0), (0, 1)
    one = monomial(0, 0)
    s_plus_one = monomial(0, 1) + one

    first = classical_k(slope) + one + one
    first_factorization = monomial(0, -1) * s_plus_one * s_plus_one
    assert first == first_factorization

    second = classical_k(add_vectors(rho, slope)) + classical_k(rho)
    second_factorization = s_plus_one * (monomial(1, 0) + monomial(-1, -1))
    assert second == second_factorization

    # Evaluation at the filled eigenvalue s=-1 kills both generators.
    for element in (first, second):
        by_r: Dict[int, int] = {}
        for (r_power, s_power), coefficient in element.terms.items():
            by_r[r_power] = by_r.get(r_power, 0) + coefficient * ((-1) ** s_power)
        assert all(value == 0 for value in by_r.values())

    return first, second


def main() -> None:
    products = exact_product_to_sum_checks()
    matrices = exact_sl2_checks()
    first, second = exact_classical_factorization_checks()

    print(f"product-to-sum identities checked: {products}")
    print(f"SL(2,Z) covariance checks: {matrices}")
    print("classical meridian relation: K_s + 2 = s^(-1)(s+1)^2")
    print("classical framing relation: K_(rho+s) + K_rho = (s+1)(rho+rho^(-1)s^(-1))")
    print(f"first classical support size: {len(first.terms)}")
    print(f"second classical support size: {len(second.terms)}")
    print("all exact checks passed")


if __name__ == "__main__":
    main()