#!/usr/bin/env python3
"""Exact integer cross-check for the rational oblique radius-six construction."""

from fractions import Fraction
from math import isqrt

DEN = 305
CROSS = 272
R2_NUM = 36 * DEN


def norm_num(m: int, n: int) -> int:
    return DEN * m * m + DEN * n * n + CROSS * m * n


def offsets(bound: int = 20) -> list[tuple[int, int]]:
    return [
        (m, n)
        for m in range(-bound, bound + 1)
        for n in range(-bound, bound + 1)
        if (m, n) != (0, 0) and norm_num(m, n) <= R2_NUM
    ]


def tri_norm(m: int, n: int) -> int:
    return m * m + m * n + n * n


obs = offsets()
assert all(abs(m) <= 6 and abs(n) <= 6 for m, n in obs)
assert offsets(6) == obs
assert min(norm_num(m, n) for m in range(-20, 21) for n in range(-20, 21) if (m, n) != (0, 0)) == 305

tri = [
    (m, n)
    for m in range(-20, 21)
    for n in range(-20, 21)
    if (m, n) != (0, 0) and tri_norm(m, n) <= 36
]
assert all(abs(m) <= 6 and abs(n) <= 6 for m, n in tri)

sum_l1 = sum(abs(m) + abs(n) for m, n in obs)
sum_prod = sum(abs(m) * abs(n) for m, n in obs)

assert len(obs) == 128
assert len(tri) == 126
assert sum_l1 == 732
assert sum_prod == 998


def directed(L: int) -> int:
    return sum((L - abs(m)) * (L - abs(n)) for m, n in obs)


assert directed(365) == 16_786_618
assert 126 * 365**2 == 16_786_350
assert directed(365) - 126 * 365**2 == 268
assert directed(365) % 2 == 0
assert directed(365) // 2 == 8_393_309
assert 63 * 365**2 == 8_393_175
assert directed(365) // 2 - 63 * 365**2 == 134
assert all(directed(L) <= 126 * L * L for L in range(7, 365))
assert directed(365) > 126 * 365**2

# Exact positive root of 2 L^2 - 732 L + 998 is between 364 and 365.
disc = 732**2 - 8 * 998
assert isqrt(disc) ** 2 <= disc < (isqrt(disc) + 1) ** 2
root = Fraction(732 + isqrt(disc), 4)  # lower rational approximation only

print(f"basis v=(136/{DEN},273/{DEN}); 136^2+273^2={136**2+273**2}={DEN**2}")
print(f"oblique offsets: {len(obs)}")
print(f"triangular offsets: {len(tri)}")
print(f"sum |m|+|n|: {sum_l1}")
print(f"sum |m||n|: {sum_prod}")
print(f"L=365 directed: {directed(365)}")
print(f"triangular directed allowance: {126*365**2}")
print(f"directed excess: {directed(365)-126*365**2}")
print(f"unordered: {directed(365)//2}")
print(f"unordered allowance: {63*365**2}")
print(f"unordered excess: {directed(365)//2-63*365**2}")
print("first winning integer L: 365")
