check_metric_group.py

#!/usr/bin/env python3
"""Exact finite regression audit for the metric-group publication.

The program checks finite rational instances of the identities used by the
proof. It deliberately does not claim to prove any infinite theorem.
"""

from __future__ import annotations

from fractions import Fraction
from itertools import product


def h(c: Fraction, x: Fraction) -> Fraction:
    return x / (c + x)


def h_inv(c: Fraction, x: Fraction) -> Fraction:
    if not 0 <= x < 1:
        raise ValueError(f"inverse input outside [0,1): {x}")
    return c * x / (1 - x) if x else Fraction(0)


def encode(c: Fraction, digits: tuple[int, ...] | list[int]) -> Fraction:
    value = Fraction(0)
    for digit in reversed(digits):
        value = Fraction(digit) + h(c, value)
    return value


def decode_finite(c: Fraction, value: Fraction, limit: int = 128) -> list[int]:
    digits: list[int] = []
    residual = value
    for _ in range(limit):
        digit = residual.numerator // residual.denominator
        digits.append(digit)
        fractional = residual - digit
        if fractional == 0:
            return digits
        residual = h_inv(c, fractional)
    raise AssertionError(f"decoder did not terminate: c={c}, value={value}")


def omega(c: Fraction, n: int) -> Fraction:
    value = Fraction(1)
    for _ in range(n):
        value = h(c, value)
    return value


def pad(values: list[int], length: int) -> list[int]:
    return values + [0] * (length - len(values))


def star_rational(c: Fraction, x: Fraction, y: Fraction) -> Fraction:
    dx = decode_finite(c, x)
    dy = decode_finite(c, y)
    length = max(len(dx), len(dy))
    return encode(c, tuple(a ^ b for a, b in zip(pad(dx, length), pad(dy, length))))


def expected_left_limit(a: list[int], q: list[int]) -> Fraction:
    k = max(i for i, digit in enumerate(a) if digit)
    prefix = [a[i] ^ q[i] for i in range(k)]
    boundary = ((a[k] - 1) ^ q[k]) + 1
    return encode(Fraction(2), tuple(prefix + [boundary]))


def observed_left_probe(a: list[int], q: list[int], large: int) -> Fraction:
    k = max(i for i, digit in enumerate(a) if digit)
    left = list(a)
    left[k] -= 1
    left[k + 1] = large
    for i in range(k + 2, len(left)):
        left[i] = 0
    return encode(Fraction(2), tuple(x ^ y for x, y in zip(left, q)))


def main() -> None:
    iterate_checks = 0
    for c in (Fraction(1), Fraction(2), Fraction(3)):
        for n in range(9):
            got = omega(c, n)
            if c == 1:
                want = Fraction(1, n + 1)
            else:
                want = (c - 1) / (c ** (n + 1) - 1)
            assert got == want, (c, n, got, want)
            iterate_checks += 1

    subadditivity_checks = 0
    words = list(product(range(4), repeat=3))
    for c in (Fraction(1, 2), Fraction(1), Fraction(2), Fraction(3)):
        for left in words:
            for right in words:
                joined = tuple(a ^ b for a, b in zip(left, right))
                assert encode(c, joined) <= encode(c, left) + encode(c, right)
                subadditivity_checks += 1

    rational_roundtrips = 0
    for c in (Fraction(1), Fraction(2), Fraction(3)):
        for denominator in range(1, 25):
            for numerator in range(3 * denominator + 1):
                value = Fraction(numerator, denominator)
                digits = decode_finite(c, value)
                assert encode(c, tuple(digits)) == value, (c, value, digits)
                rational_roundtrips += 1

    sample_identities = {
        "half_star_one": star_rational(Fraction(2), Fraction(1, 2), Fraction(1)),
        "one_star_four_thirds": star_rational(
            Fraction(2), Fraction(1), Fraction(4, 3)
        ),
        "star1_separator": star_rational(
            Fraction(1), Fraction(2, 5), Fraction(1, 3)
        ),
        "star2_separator": star_rational(
            Fraction(2), Fraction(2, 5), Fraction(1, 3)
        ),
    }
    assert sample_identities == {
        "half_star_one": Fraction(3, 2),
        "one_star_four_thirds": Fraction(1, 3),
        "star1_separator": Fraction(3, 7),
        "star2_separator": Fraction(1, 7),
    }

    c_half = Fraction(1, 2)
    gap_checks = 0
    for position in range(12):
        label = encode(c_half, tuple([0] * position + [1]))
        assert label > Fraction(1, 2)
        gap_checks += 1
    assert h(c_half, Fraction(1, 2)) == Fraction(1, 2)

    boundary_probes = 0
    boundary_failures = 0
    finite_points = [
        [a0, a1, a2, 0, 0]
        for a0 in range(5)
        for a1 in range(5)
        for a2 in range(5)
        if any((a0, a1, a2))
    ]
    translators = [
        [q0, q1, q2, 0, 0]
        for q0 in range(5)
        for q1 in range(5)
        for q2 in range(5)
    ]
    for point in finite_points:
        k = max(i for i, digit in enumerate(point) if digit)
        for translator in translators:
            boundary_probes += 1
            predicted = expected_left_limit(point, translator)
            observed = observed_left_probe(point, translator, 10**8)
            if abs(observed - predicted) >= Fraction(1, 10**6):
                boundary_failures += 1
            direct = encode(
                Fraction(2), tuple(x ^ y for x, y in zip(point, translator))
            )
            tail_zero = all(digit == 0 for digit in translator[k + 1 :])
            m = point[k]
            cancellation = tail_zero and (
                (m ^ translator[k]) == (((m - 1) ^ translator[k]) + 1)
            )
            if cancellation != (predicted == direct):
                boundary_failures += 1
    assert boundary_failures == 0

    print(f"iterate identities: {iterate_checks}")
    print(f"finite XOR-subadditivity pairs: {subadditivity_checks}")
    print(f"integer-parameter rational roundtrips: {rational_roundtrips}")
    print(f"sharp-threshold gap probes: {gap_checks}")
    print(f"translation-boundary probes: {boundary_probes}")
    print(f"translation-boundary failures: {boundary_failures}")
    print("metric-group finite audit: PASS")


if __name__ == "__main__":
    main()