#!/usr/bin/env python3
"""Exact synthesis certificate for the all-rank degree-six theorem."""

from __future__ import annotations

import hashlib
import subprocess
import sys
from fractions import Fraction
from itertools import combinations
from math import factorial
from pathlib import Path

import sympy as sp


def repository_root(path: Path) -> Path:
    for parent in path.parents:
        if (parent / "bin" / "rhc.py").is_file():
            return parent
    raise RuntimeError("repository root not found")


ROOT = repository_root(Path(__file__).resolve())
DEPENDENCIES = {
    "canon/claims/C-0043-sharp-degree-thirteen-karlin-wronskian.md":
        "060a5c3956ce337ec0d065c1aa3a52d11c8358529b596588a30f6f4937d7ba06",
    "canon/witnesses/C-0043/verify.py":
        "5049666042b70c35fec02bddf6a22deb1ea797239ae7f8bfffc43ffcce915e09",
    "canon/claims/C-0056-sharp-degree-eight-karlin-wronskian.md":
        "e5c24cab55e3bf6ce6e2caf9b10542cbbf61ce6986b3e83c9d1a4d45c70d41af",
    "canon/witnesses/C-0056/verify.py":
        "81705c43d0db60f5d5dbed189ae8220bb4943cf307407cf74f3bee82c8ec4623",
    "canon/claims/C-0059-sharp-degree-seven-signed-karlin-wronskian.md":
        "8e540fa522da32b346f616bfc2a66c5d873ac716638d7d80801a2c0de72333ec",
    "canon/witnesses/C-0059/verify.py":
        "5406185dff46e2aae06a5a55ecd3d8afb4dd4a4c9a204cbc74bc4f44f29115d3",
}
DEPENDENCY_VERIFIERS = (
    "canon/witnesses/C-0043/verify.py",
    "canon/witnesses/C-0056/verify.py",
    "canon/witnesses/C-0059/verify.py",
)
EXPECTED_BOUNDARY_WITNESS = -1_246_946_918_400
EXPECTED_SIMPLE_WITNESS_SHA256 = (
    "3cee62a7790b6341ed2b586e6cf3dd0add9e3e9499bca2116bc1846f2569c5ca"
)


if not __debug__:
    raise SystemExit("verify.py requires assertions")


def sha256(path: Path) -> str:
    return hashlib.sha256(path.read_bytes()).hexdigest()


def verify_dependencies() -> None:
    for relative, expected in DEPENDENCIES.items():
        assert sha256(ROOT / relative) == expected

    for relative in DEPENDENCY_VERIFIERS:
        completed = subprocess.run(
            [sys.executable, str(ROOT / relative)],
            cwd=ROOT,
            check=False,
            capture_output=True,
            text=True,
        )
        assert completed.returncode == 0, completed.stderr
        assert "PASS" in completed.stdout


def elementary_symmetric(
    values: tuple[sp.Expr | Fraction, ...],
    degree: int,
) -> sp.Expr:
    return sp.Add(
        *(
            sp.prod(values[index] for index in subset)
            for subset in combinations(range(len(values)), degree)
        )
    )


def product_derivatives(
    values: tuple[sp.Expr | Fraction, ...],
    maximum: int,
) -> list[sp.Expr]:
    return [
        factorial(degree) * elementary_symmetric(values, degree)
        if degree <= len(values)
        else sp.Integer(0)
        for degree in range(maximum + 1)
    ]


def signed_hankel(derivatives: list[sp.Expr], rank: int) -> sp.Expr:
    determinant = sp.Matrix(
        [
            [
                derivatives[row + column]
                for column in range(rank)
            ]
            for row in range(rank)
        ]
    ).det()
    return sp.expand(
        (-1) ** (rank * (rank - 1) // 2) * determinant
    )


def verify_low_ranks() -> None:
    alpha = sp.symbols("a0:6")
    derivatives = product_derivatives(alpha, 4)
    power = {
        degree: sum(value**degree for value in alpha)
        for degree in range(1, 5)
    }

    assert signed_hankel(derivatives, 1) == 1
    assert sp.expand(
        signed_hankel(derivatives, 2) - power[2]
    ) == 0

    rank_three = 2 * (
        power[2] ** 3
        - 3 * power[2] * power[4]
        + 2 * power[3] ** 2
    )
    assert sp.expand(
        signed_hankel(derivatives, 3) - rank_three
    ) == 0

    positive_formula = 8 * sum(
        alpha[left] ** 3 * alpha[right] ** 3
        for left, right in combinations(range(6), 2)
    ) + 12 * sum(
        sp.prod(alpha[index] ** 2 for index in subset)
        for subset in combinations(range(6), 3)
    )
    assert sp.expand(rank_three - positive_formula) == 0


def verify_high_ranks() -> None:
    coefficients = sp.symbols("c0:7")
    derivatives = [
        factorial(degree) * coefficients[degree]
        if degree <= 6
        else sp.Integer(0)
        for degree in range(19)
    ]
    expected_rank_seven = (factorial(6) * coefficients[6]) ** 7
    assert sp.expand(
        signed_hankel(derivatives, 7) - expected_rank_seven
    ) == 0
    for rank in range(8, 11):
        assert signed_hankel(derivatives, rank) == 0


def centered_determinant(
    values: tuple[Fraction, ...],
    rank: int,
) -> sp.Expr:
    derivatives = product_derivatives(values, 2 * rank - 2)
    return signed_hankel(derivatives, rank)


def verify_degree_seven_counterfeit() -> None:
    # f(x)=x^4(x+1)^3, represented by ascending coefficients.
    coefficients = (0, 0, 0, 0, 1, 3, 3, 1)
    derivatives = [
        factorial(degree) * (
            coefficients[degree] if degree < len(coefficients) else 0
        )
        for degree in range(11)
    ]
    assert signed_hankel(derivatives, 6) == EXPECTED_BOUNDARY_WITNESS

    denominator = 10_000_000
    numerators = (
        53 * denominator - 3,
        53 * denominator - 1,
        53 * denominator + 1,
        53 * denominator + 3,
        denominator - 2,
        denominator,
        denominator + 2,
    )
    reciprocal_roots = tuple(
        Fraction(numerator, denominator)
        for numerator in numerators
    )
    assert len(set(reciprocal_roots)) == 7
    assert min(reciprocal_roots) > 0
    determinant = Fraction(centered_determinant(reciprocal_roots, 6))
    payload = f"{determinant.numerator}/{determinant.denominator}".encode()
    assert hashlib.sha256(payload).hexdigest() == (
        EXPECTED_SIMPLE_WITNESS_SHA256
    )
    assert determinant < 0


def main() -> None:
    verify_dependencies()
    verify_low_ranks()
    verify_high_ranks()
    verify_degree_seven_counterfeit()
    print("C-0043, C-0056, and C-0059 hashes and verifiers: PASS")
    print("ranks one through three exact formulas: PASS")
    print("rank-seven anti-triangular formula: PASS")
    print("ranks eight and higher zero-column reduction: PASS")
    print("degree-seven boundary and simple-root counterfeits: PASS")


if __name__ == "__main__":
    main()
