2607.19268v1 / tools/check_finite.py

all files

#!/usr/bin/env python3
"""Independent numerical audit of plantri's finite output.

This is deliberately not part of the Lean proof.  It reads graph6 lines
produced by plantri, applies power iteration, and reports the two largest
adjacency spectral radii together with the value for K₂ ∨ P_(n-2).
"""

from __future__ import annotations

import argparse
import math
from pathlib import Path


def decode_graph6(line: bytes) -> list[list[int]]:
    line = line.strip()
    if line.startswith(b">>graph6<<"):
        line = line[len(b">>graph6<<") :]
    if not line or line[0] == 126:
        raise ValueError("the audit only implements the one-byte graph6 order")
    n = line[0] - 63
    bits: list[int] = []
    for byte in line[1:]:
        value = byte - 63
        bits.extend((value >> shift) & 1 for shift in range(5, -1, -1))
    adjacency = [[] for _ in range(n)]
    cursor = 0
    for upper in range(1, n):
        for lower in range(upper):
            if bits[cursor]:
                adjacency[lower].append(upper)
                adjacency[upper].append(lower)
            cursor += 1
    return adjacency


def join_k2_path(n: int) -> list[list[int]]:
    if n < 3:
        raise ValueError("n must be at least 3")
    edges = {(0, 1)}
    for vertex in range(2, n):
        edges.add((0, vertex))
        edges.add((1, vertex))
    for vertex in range(2, n - 1):
        edges.add((vertex, vertex + 1))
    adjacency = [[] for _ in range(n)]
    for left, right in edges:
        adjacency[left].append(right)
        adjacency[right].append(left)
    return adjacency


def spectral_radius(adjacency: list[list[int]], iterations: int = 500) -> float:
    n = len(adjacency)
    vector = [1.0 / math.sqrt(n)] * n
    quotient = 0.0
    for _ in range(iterations):
        product = [sum(vector[j] for j in adjacency[i]) for i in range(n)]
        norm = math.sqrt(sum(value * value for value in product))
        vector = [value / norm for value in product]
        next_product = [sum(vector[j] for j in adjacency[i]) for i in range(n)]
        next_quotient = sum(
            vector[i] * next_product[i] for i in range(n)
        )
        if abs(next_quotient - quotient) < 1e-14:
            return next_quotient
        quotient = next_quotient
    return quotient


def main() -> None:
    parser = argparse.ArgumentParser()
    parser.add_argument("n", type=int)
    parser.add_argument("graph6", type=Path)
    args = parser.parse_args()

    best: list[tuple[float, bytes, list[int]]] = []
    count = 0
    with args.graph6.open("rb") as stream:
        for raw in stream:
            if raw.startswith(b">>"):
                continue
            adjacency = decode_graph6(raw)
            if len(adjacency) != args.n:
                raise ValueError(f"expected order {args.n}, got {len(adjacency)}")
            radius = spectral_radius(adjacency)
            degrees = sorted((len(row) for row in adjacency), reverse=True)
            best.append((radius, raw.strip(), degrees))
            best.sort(reverse=True)
            del best[3:]
            count += 1

    target = spectral_radius(join_k2_path(args.n))
    print(f"n={args.n} count={count} target={target:.12f}")
    for rank, (radius, code, degrees) in enumerate(best, 1):
        print(
            f"  {rank}: radius={radius:.12f} "
            f"delta={radius - target:+.12e} degrees={degrees} graph6={code!r}"
        )


if __name__ == "__main__":
    main()