"""Distance spectral radii of connected graphs -- numberdb.org/T536.

Run with SageMath and the attached files in one directory:
    sage -pip install numberdb
    sage -python generate.py             # verify all stored rows (read only)
    sage -python checks.py               # offline independent build checks

To fill this draft with NUMBERDB_API_KEY configured:
    sage -python generate.py --fill
This writes entries, never offers, reviews, or publishes the website draft.

The graph6 identities are frozen in selected-graphs.json, matching the Sage
canonical_label(algorithm='sage') keys of T237 and T238. They contain every
connected isomorphism class of orders 2 through 6. checks.py independently
compares this set with nauty geng and McKay's catalogues.

All matrices and characteristic polynomials are exact over ZZ. AA isolates
all real roots with multiplicity and selects the largest exactly. A rational
root is returned exactly. Other roots are enclosed in a real interval and
then a real ball at bits(digits, losing=64). The measured width and the
independent rational Collatz--Wielandt checks are in numerical-evidence.json.
"""
import json
import os
from pathlib import Path
import re
import sys
from functools import lru_cache

import numberdb.sage as numberdb
from sage.graphs.graph import Graph
from sage.matrix.constructor import matrix
from sage.rings.integer_ring import ZZ
from sage.rings.rational_field import QQ
from sage.rings.qqbar import AA
from sage.rings.real_mpfi import RealIntervalField
from sage.rings.real_arb import RealBallField

WORKING_GUARD = 64
HERE = Path(__file__).resolve().parent


def selected_graphs():
    rows = json.loads((HERE / 'selected-graphs.json').read_text())
    keys = [row['g'] for row in rows]
    assert len(keys) == len(set(keys)) == 142
    assert keys == sorted(keys, key=lambda g: (ord(g[0]) - 63, g))
    assert {n: sum(ord(g[0])-63 == n for g in keys) for n in range(2, 7)} == {
        2: 1, 3: 2, 4: 6, 5: 21, 6: 112}
    return rows


@lru_cache(maxsize=None)
def spectral_data(key):
    graph = Graph(key, format='graph6')
    assert graph.is_connected() and graph.order() >= 2
    distance = matrix(ZZ, graph.distance_matrix())
    polynomial = distance.charpoly('x')
    roots = polynomial.roots(ring=AA)
    assert sum(m for _, m in roots) == graph.order()
    largest = max(root for root, _ in roots)
    assert largest > 0
    return distance, polynomial, roots, largest


class DistanceSpectralRadii(numberdb.Generator):
    table = 'T536'
    parameters = ('g',)
    type = 'R'
    digits = 100
    rigour = 'proven'
    files = ('generate.py', 'checks.py', 'check.py', 'selected-graphs.json',
             'argument-selection.md', 'numerical-evidence.json',
             'skill-evidence.json', 'source-notes.md', 'source-files.json',
             'fiedler-reference.json', 'graph2c.g6', 'graph3c.g6', 'graph4c.g6', 'graph5c.g6', 'graph6c.g6')

    def enumerate(self):
        for row in selected_graphs():
            yield {'g': row['g']}

    def value(self, params, digits):
        key = params['g']
        _, _, _, largest = spectral_data(key)
        if largest in QQ:
            value = QQ(largest)
        else:
            precision = numberdb.bits(digits, losing=WORKING_GUARD)
            interval = RealIntervalField(precision)(largest)
            value = RealBallField(precision)(interval)
            assert value.is_finite()
            assert value.parent()(value.rad()) < value.abs() * QQ(1) / ZZ(10)**(digits + 10)
        entry = {'number': value}
        comment = next((r.get('comment') for r in selected_graphs() if r['g'] == key), None)
        if comment:
            entry['comment'] = comment
        return entry


if __name__ == '__main__':
    if os.environ.get('NUMBERDB_KEY_FROM_STDIN') == '1':
        key = sys.stdin.read().strip()
        if re.match(r'^[A-Z_][A-Z0-9_]*\s*=', key):
            key = key.split('=', 1)[1].strip().strip('\"\'')
        numberdb.configure(api_key=key)
    generator = DistanceSpectralRadii()
    if '--fill' in sys.argv or os.environ.get('NUMBERDB_PUBLISH') == '1':
        print(generator.publish(message='Fill 142 distance radii with certified roots and independent rational checks'))
    report = generator.verify(sample=None)
    print(report)
    sys.exit(0 if report.ok else 1)
