#!/usr/bin/env python3
"""Verify the packaged exact counts and rebuild the release report."""
from __future__ import annotations

import argparse
from decimal import Decimal, getcontext
from pathlib import Path


def load(path: Path) -> list[int]:
    values: list[int] = []
    for raw in path.read_text().splitlines():
        line = raw.strip()
        if not line or line.startswith('#'):
            continue
        k_s, value_s = line.split()
        k = int(k_s)
        if k != len(values):
            raise ValueError(f'unexpected index {k} in {path}')
        values.append(int(value_s))
    return values


def partition_numbers(nmax: int) -> list[int]:
    p = [0] * (nmax + 1)
    p[0] = 1
    for n in range(1, nmax + 1):
        total = 0
        k = 1
        while True:
            g1 = k * (3 * k - 1) // 2
            if g1 > n:
                break
            sign = 1 if k % 2 else -1
            total += sign * p[n - g1]
            g2 = k * (3 * k + 1) // 2
            if g2 <= n:
                total += sign * p[n - g2]
            k += 1
        p[n] = total
    return p


def main() -> None:
    parser = argparse.ArgumentParser()
    parser.add_argument('--counts', type=Path, required=True)
    parser.add_argument('--public-snapshot', type=Path, required=True)
    parser.add_argument('--output', type=Path, required=True)
    args = parser.parse_args()

    counts = load(args.counts)
    public = load(args.public_snapshot)
    if len(counts) != 10001:
        raise AssertionError(f'expected 10001 counts, found {len(counts)}')
    if len(public) != 751:
        raise AssertionError(f'expected 751 public terms, found {len(public)}')
    if counts[:751] != public:
        raise AssertionError('packaged counts do not match the public b-file snapshot')

    nmax = 2 * (len(counts) - 1)
    p = partition_numbers(nmax)
    for k, value in enumerate(counts):
        if not 0 <= value <= p[2 * k]:
            raise AssertionError(f'count outside [0,p(2k)] at k={k}')

    getcontext().prec = 40
    density = Decimal(counts[-1]) / Decimal(p[nmax])
    report = (
        'Exact count verification\n'
        '========================\n'
        f'Range: k=0,...,{len(counts)-1}, corresponding to even total sizes N=0,...,{nmax}.\n'
        f'Number of reconstructed values: {len(counts)}.\n'
        'The exact run used 17 prime moduli with a 527-bit CRT product.\n'
        f'p({nmax}) has {p[nmax].bit_length()} bits.\n'
        'Since 0 <= a_N <= p(N), the CRT reconstruction is unique.\n'
        'The original run checked every reconstructed value modulo the independent prime 2147483171.\n'
        'The packaged values for k=0,...,750 agree with the preceding public OEIS b-file snapshot.\n'
        f'a_{nmax} = {counts[-1]}\n'
        f'a_{nmax} / p({nmax}) = {density:.24E}\n'
    )
    args.output.write_text(report)
    print(args.output)


if __name__ == '__main__':
    main()
