"""Exact checks for the prose repair of numberdb.org/T498.

Run with SageMath after downloading generate.py from the table:
    sage -python repair_checks.py
Pass an optional API JSON document filename to check its stored intervals too.
No writes are made. Catalogue identities for Lambda_1,...,Lambda_8 are sourced
from https://www.math.rwth-aachen.de/~Gabriele.Nebe/LATTICES/ ; their named
representatives are A1,A2,A3,D4,D5,E6,E7,E8 up to similarity.
"""
import json
import sys
from decimal import Decimal
from pathlib import Path
import numberdb.sage as numberdb
from sage.rings.rational_field import QQ
from sage.matrix.constructor import matrix
from sage.rings.polynomial.polynomial_ring_constructor import PolynomialRing
from sage.env import SAGE_VERSION
from generate import exact_form


def rows(m):
    return [[str(x) for x in row] for row in m.rows()]


def rational(text):
    if '/' in text:
        a, b = text.split('/')
        return QQ(int(a)) / QQ(int(b))
    a, b = Decimal(text).as_integer_ratio()
    return QQ(a) / QQ(b)


def main():
    evidence = {'sage_version': SAGE_VERSION, 'similarities': {}}
    a1 = matrix(QQ, [[2]])
    assert a1.inverse() == a1 / 4
    evidence['similarities']['A1*'] = 'Gram(A1*) = Gram(A1)/4'
    a2 = matrix(QQ, [[2,-1],[-1,2]])
    u = matrix(QQ, [[1,0],[0,-1]])
    assert abs(u.det()) == 1 and a2.inverse() == u*a2*u.transpose()/3
    evidence['similarities']['A2*'] = 'Gram(A2*) = diag(1,-1) Gram(A2) diag(1,-1)/3'
    d4 = matrix(QQ, [[1,-1,0,0],[0,1,-1,0],[0,0,1,-1],[0,0,1,1]])
    transform = matrix(QQ, [[1,1,0,0],[-1,1,0,0],[0,0,1,1],[0,0,-1,1]])/2
    eye4 = matrix(QQ, 4, 4, lambda i,j: int(i==j))
    assert transform*transform.transpose() == eye4/2
    change = d4*transform*d4.transpose()
    assert all(x.denominator()==1 for x in change.list()) and abs(change.det())==1
    assert change*d4.inverse().transpose() == d4*transform
    evidence['similarities']['D4*'] = {'squared_scale':'1/2','unimodular_change':rows(change)}
    e8 = matrix(QQ, [[2,0,0,0,0,0,0,0],[-1,1,0,0,0,0,0,0],
        [0,-1,1,0,0,0,0,0],[0,0,-1,1,0,0,0,0],
        [0,0,0,-1,1,0,0,0],[0,0,0,0,-1,1,0,0],
        [0,0,0,0,0,-1,1,0],[QQ(1)/2]*8])
    gram=e8*e8.transpose()
    assert gram.det()==1 and all(x.denominator()==1 for x in gram.list())
    assert all(x.denominator()==1 for x in gram.inverse().list())
    assert gram.inverse()*e8 == e8.inverse().transpose()
    evidence['similarities']['E8*'] = 'Integral Gram matrix of determinant 1 proves E8 is self-dual.'
    d3=matrix(QQ, [[1,-1,0],[0,1,-1],[0,1,1]])
    d3=d3.matrix_from_rows([1,0,2])
    a3=matrix(QQ, [[2,-1,0],[-1,2,-1],[0,-1,2]])
    assert d3*d3.transpose()==a3
    evidence['similarities']['D3=A3'] = 'Equal Gram matrices after permuting the D3 basis.'
    for f,n in [('A',1),('A',2),('D',4)]:
        a,b=exact_form(f,n); ad,bd=exact_form(f+'*',n)
        assert a**n*b == ad**n*bd
    evidence['omitted_dual_moments'] = 'Exact positive nth powers agree for A1/A1*, A2/A2*, D4/D4*.'
    ring=PolynomialRing(QQ, names=('t','s'))
    t,s=ring.gens(); field=ring.fraction_field(); t=field(t); s=field(s)
    for n in range(1,25):
        assert s/(n*t**(n+2)) == (s/t**n)/(n*t**2)
    evidence['conversion'] = 'With v=t^n, I=s/t^n and G=s/(n*t^(n+2)); conversion checked exactly for n=1..24.'
    # The general identity follows from v^(1+2/n)=v*v^(2/n), for v>0.
    if len(sys.argv)>1:
        path=Path(sys.argv[1])
    else:
        path=Path(__file__).with_name('T498-before.json')
    if path.exists():
        doc=json.loads(path.read_text())
        count=0; independent=0
        for family, entries in doc['Numbers'].items():
            for ns,item in entries.items():
                n=int(ns); text=item['number'] if isinstance(item,dict) else item
                a,b=exact_form(family,n); expected=a**n*b
                center=rational(text)
                if '.' in text or 'e' in text.lower():
                    ulp=QQ(10)**Decimal(text).as_tuple().exponent
                    assert center-ulp>0
                    assert (center-ulp)**n <= expected <= (center+ulp)**n, (family,n)
                else:
                    assert center**n == expected, (family,n)
                if family in ('Z','A','D','E'):
                    # Independent transcription of CS82 equations (32)-(36).
                    if family=='Z': target=(QQ(1)/12)**n
                    elif family=='A': target=((QQ(n+3)/(12*(n+1)))**n)/(n+1)
                    elif family=='D': target=((QQ(n*(n+1)+6)/(12*n*(n+1)))**n)/4
                    else: target={6:(QQ(5)/56)**6/3,7:(QQ(163)/2016)**7/2,8:(QQ(929)/12960)**8}[n]
                    assert expected==target
                    independent+=1
                count+=1
        assert count==114 and independent==72
        evidence['stored_entries_checked_exactly']=count
        evidence['entries_also_compared_with_source_formulas']=independent
    print(json.dumps(evidence,indent=2))


if __name__=='__main__':
    main()
