"""Covering radii of spherical coverings -- numberdb.org/T480.

Reproduce with SageMath and the attached source.json:
    sage -pip install numberdb
    sage -python generate.py                 # verify every stored entry
    sage -python checks.py                   # independent geometry checks
    sage -python generate.py --fill-draft    # fill existing draft, key required

Set NUMBERDB_API_KEY for authenticated reads of the draft. This script does
not create tables, offer them for review, or make them public. The optional
NUMBERDB_KEY_FROM_STDIN=1 reads the key from stdin instead.

Hardin, Sloane and Smith's 1994 configurations are taken from Henry Cohn's
DSpace archive, downloaded 2026-10-25. The PDF rounds angular radii UP to
seven decimal places. The source model is [printed - 1e-7, printed] degrees.
This numerical uncertainty is not a bound on unknown global optima.

Seven elementary configurations use their exact geometry at 100 significant
digits. The remaining angles retain the source's seven decimal places;
chords use eight significant digits, checked to enclose the transformed
source interval. WORKING_GUARD=80 leaves over 100 decimal digits of accuracy
in every elementary result, as measured by checks.py at two precisions.

The polar configurations for n=10,14 consist of two poles and two staggered
regular m-gons, m=(n-2)/2, at heights +/-c/(1+c), c=cos(pi/m). Their equal
triangular facet distances give the formula implemented here. The attached
checks construct these points separately and verify every hull facet.
"""
import numberdb.sage as numberdb
from sage.rings.real_arb import RealBallField
from sage.rings.rational_field import QQ
from decimal import Decimal, localcontext
import json
import os
import pathlib
import sys

WORKING_GUARD = 80
EXACT = frozenset((4, 5, 6, 7, 10, 12, 14))
SOURCE = json.loads(pathlib.Path(__file__).with_name('source.json').read_text())


def rational_decimal(text):
    """Parse the finite decimal as a rational, without binary floating point."""
    whole, fraction = text.split('.')
    return QQ(int(whole + fraction)) / QQ(10 ** len(fraction))


def cosine(n, field):
    if n == 4:
        return field(1) / 3
    if n in (5, 6):
        return 1 / field(5 if n == 5 else 3).sqrt()
    if n == 7:
        c = (1 + field(5).sqrt()) / 4
        return c / (1 + c*c).sqrt()
    if n == 12:
        return ((5 + 2*field(5).sqrt()) / 15).sqrt()
    if n in (10, 14):
        c = field(2 if n == 10 else 3).sqrt() / 2
        t = c*c*(1 + 2*c)
        return (t / (1+t)).sqrt()
    raise ValueError('No elementary formula specified for this n.')


def source_angle(n, field):
    printed = rational_decimal(SOURCE['rows'][str(n)]['angle'])
    half_ulp = QQ(1) / (2 * 10**7)
    return field(printed-half_ulp, half_ulp)


def entry_comment(n):
    formulas = {
        4: r'Regular tetrahedron: $\cos(\pi\rho_4/180)=1/3$. Optimality is proven CITE{Fejes1943}.',
        5: r'Triangular bipyramid: $\cos(\pi\rho_5/180)=1/\sqrt5$. Optimality is proven CITE{Schutte}.',
        6: r'Regular octahedron: $\cos(\pi\rho_6/180)=1/\sqrt3$. Optimality is proven CITE{Fejes1943}.',
        7: r'Pentagonal bipyramid: $\cos(\pi\rho_7/180)=\sqrt{(3+\sqrt5)/(11+\sqrt5)}$. Optimality is proven CITE{Schutte}.',
        8: r"Schütte's arrangement. Optimality is proven CITE{Wimmer}; these entries retain numerical source precision.",
        10: r'Two poles and two staggered squares: $\cos(\pi\rho_{10}/180)=\sqrt{(1+\sqrt2)/(3+\sqrt2)}$. Optimality is proven CITE{Fejes1969}.',
        12: r'Regular icosahedron: $\cos(\pi\rho_{12}/180)=\sqrt{(5+2\sqrt5)/15}$. Optimality is proven CITE{Fejes1943}.',
        14: r'Two poles and two staggered regular hexagons: $\cos(\pi\rho_{14}/180)=\sqrt{(3+3\sqrt3)/(7+3\sqrt3)}$. Optimality is proven CITE{Fejes1969}.',
    }
    return formulas.get(n)


class SphericalCoverings(numberdb.Generator):
    table = 'T480'
    parameters = ('n', 'normalisation')
    type = 'R'
    digits = 100
    rigour = 'heuristic (agreement-checked)'
    files = ('generate.py', 'checks.py', 'source.json', 'check.py',
             'evidence.json', 'README.txt')

    def enumerate(self):
        for n in range(4, 131):
            for normalisation in ('angle', 'chord'):
                yield {'n': n, 'normalisation': normalisation}

    def digits_for(self, params):
        if int(params['n']) in EXACT:
            return 100
        return 9 if params['normalisation'] == 'angle' else 8

    def value(self, params, digits):
        n = int(params['n'])
        normalisation = params['normalisation']
        field = RealBallField(numberdb.bits(max(100, digits), losing=WORKING_GUARD))
        if n in EXACT:
            c = cosine(n, field)
            value = c.arccos()*180/field.pi() if normalisation == 'angle' else (2-2*c).sqrt()
        elif normalisation == 'angle':
            value = SOURCE['rows'][str(n)]['angle']
        else:
            ball = 2*(source_angle(n, field)*field.pi()/360).sin()
            with localcontext() as context:
                context.prec = 8
                value = format(+Decimal(str(ball.mid())), 'f')
            centre = rational_decimal(value)
            ulp = QQ(1) / 10**len(value.split('.')[1])
            assert field(centre-ulp) < ball < field(centre+ulp), (n,value,str(ball),str(ulp))
        result = {'number': value, 'digits': self.digits_for(params)}
        comment = entry_comment(n)
        if comment:
            result['comment'] = comment
        return result


if __name__ == '__main__':
    if os.environ.get('NUMBERDB_KEY_FROM_STDIN') == '1':
        os.environ['NUMBERDB_API_KEY'] = sys.stdin.read().strip()
    generator = SphericalCoverings()
    if '--fill-draft' in sys.argv or os.environ.get('NUMBERDB_PUBLISH') == '1':
        outcome = generator.publish(message='Fill spherical covering radii from archived coordinates and exact small configurations.')
        print(outcome)
        print('FILL_RUN', outcome.run)
        print('FILES', outcome.files)
    else:
        report = generator.verify(sample=None)
        print(report)
        sys.exit(0 if report.ok else 1)