"""Independent checks for T662; run: sage -python checks.py [--stored].

MPFI checks signs at both endpoints of every displayed decimal interval using
the original cost equation. mpmath secant solves give a second numerical path.
Direct finite group multiplication checks CKSU's triple-product subsets at
m=2,3,4,5; a deliberately invalid construction is a negative control.
GAP character tables check the degree multiplicities at m=2,3 after an S3
control. The live companion T661 supplies all 36 cubic and quarter-power costs.
None of these checks writes to the website. --stored repeats the numerical
checks on the actual stored text after filling. Numerical evidence identifies
the archived skill manifest; this is build evidence, not independent critique.
"""
import sys,json,platform
from fractions import Fraction
from itertools import product
from collections import Counter
import numberdb.sage as numberdb
from sage.rings.integer_ring import ZZ
from sage.rings.rational_field import QQ
from sage.rings.real_mpfi import RealIntervalField
from sage.rings.real_arb import RealBallField
from sage.env import SAGE_VERSION
from numberdb._write import to_text
import mpmath as mp
from generate import ConstructionBounds, certified_bisect


def endpoints(text):
    center = Fraction(text)
    digits_after = len(text.split('.')[1])
    ulp = Fraction(1,10**digits_after)
    return QQ((center-ulp).numerator)/QQ((center-ulp).denominator), QQ((center+ulp).numerator)/QQ((center+ulp).denominator)


def mul(x,y,m):
    a,b,e=x; c,d,f=y
    if e: c,d=d,c
    return (tuple((a[i]+c[i])%m for i in range(3)),
            tuple((b[i]+d[i])%m for i in range(3)),e^f)


def inv(x,m):
    a,b,e=x
    if e: a,b=b,a
    return tuple(-v%m for v in a),tuple(-v%m for v in b),e


def subsets(m,include_zero=False):
    sets=[]
    for i in range(3):
        subset=[]
        for a,b,e in product(range(0 if include_zero else 1,m),range(m),range(2)):
            av=[0]*3; bv=[0]*3; av[i]=a; bv[(i+1)%3]=b
            subset.append((tuple(av),tuple(bv),e))
        sets.append(subset)
    return sets


def violations(m,include_zero=False):
    sets=subsets(m,include_zero)
    quotients=[{mul(x,inv(y,m),m) for x in subset for y in subset} for subset in sets]
    identity=((0,0,0),(0,0,0),0)
    count=0
    for x in quotients[0]:
        for y in quotients[1]:
            z=inv(mul(x,y,m),m)
            if z in quotients[2] and (x!=identity or y!=identity or z!=identity): count+=1
    return count,[len(s) for s in sets],[len(q) for q in quotients]


def main():
    R=RealIntervalField(512); B=RealBallField(429)
    mp.mp.dps=150
    # Controls precede the measurements.
    lo,hi=certified_bisect(lambda x:B(x)**2-2,1,2,357)
    assert lo < R(2).sqrt().lower() and R(2).sqrt().upper() < hi
    control=mp.findroot(lambda x:x*x-2,(1,2))
    assert abs(control-mp.sqrt(2)) < mp.mpf('1e-145')
    bad,_,_=violations(2,True); assert bad>0
    print('CONTROLS',json.dumps({'sqrt2_arb_bracket':True,'sqrt2_secant':True,'invalid_TPP_violations':bad}),flush=True)
    doc=numberdb.table('T662') if '--stored' in sys.argv else None
    companion=numberdb.table('T661')
    results=[]; texts={}; min_accuracy=10000
    for params in ConstructionBounds().enumerate():
        m=int(params['m']); q=m**3; A=2*m*(m-1)
        ball=ConstructionBounds().value(params,100)
        assert ball.is_finite() and QQ(ball.rad()) < QQ(1)/10**105
        min_accuracy=min(min_accuracy,int((QQ(2)/QQ(ball.rad())).floor()).bit_length()-1)
        text=doc['Numbers'][str(m)] if doc else to_text(ball,100)
        if isinstance(text,dict): text=text['number']
        low,high=endpoints(text)
        def original(s):
            return (R(A).log()*R(s)).exp()-2*q-R(q*(q-1)//2)*(R(2).log()*R(s)).exp()
        left,right=original(low),original(high)
        assert left.upper()<0 and right.lower()>0, (m,'display bracket')
        root=mp.findroot(lambda s:mp.power(A,s)-2*q-mp.mpf(q*(q-1)//2)*mp.power(2,s),(mp.mpf('2.8'),mp.mpf('3.1')))
        mp_low=mp.mpf(str(low.numerator()))/mp.mpf(str(low.denominator()))
        mp_high=mp.mpf(str(high.numerator()))/mp.mpf(str(high.denominator()))
        assert mp_low < root < mp_high
        coarse=(R(m)*0+6*R(m).log()-R(2).log())/R(m*(m-1)).log()
        assert R(high)<coarse and coarse<3
        # Independent dimension-square identity and moments read from companion.
        assert 2*q+4*(q*(q-1)//2)==2*m**6
        cubic=companion['Numbers'][str(q)]['2']['3']
        assert int(cubic)==4*q*q-2*q
        quarter=companion['Numbers'][str(q)]['2']['9/4']
        cl,ch=endpoints(quarter)
        cost=2*q+R(q*(q-1)//2)*(R(2).log()*R(QQ(9)/4)).exp()
        assert R(cl)<cost and cost<R(ch)
        assert 2<low<high<3
        texts[str(m)]=text
        results.append({'m':m,'value':text,'display_endpoint_signs':['negative','positive'],
                        'independent_secant_inside':True,'coarse_upper_bound':str(coarse)})
    minimum=min(results,key=lambda row:Fraction(row['value']))
    assert minimum['m']==17
    print('NUMERICAL',json.dumps({'mode':'stored' if doc else 'pre-fill','entries':len(results),'min_relative_accuracy_bits':min_accuracy,'minimum_m_in_selection':minimum['m'],'maximum_display_length':max(map(len,texts.values())),'display_block_bytes':len(json.dumps(texts))}),flush=True)
    print('VALUES',json.dumps(results),flush=True)
    for m in (2,3,4,5):
        count,sizes,qsizes=violations(m)
        assert count==0 and sizes==[2*m*(m-1)]*3
        print('TPP',json.dumps({'m':m,'subset_sizes':sizes,'quotient_sizes':qsizes,'violations':count}),flush=True)
    from sage.interfaces.gap import gap
    import ast
    degrees=ast.literal_eval(gap.eval('List(Irr(SymmetricGroup(3)), chi -> chi[1])'))
    assert sorted(degrees)==[1,1,2]
    for m in (2,3):
        command='Collected(List(Irr(WreathProduct(AbelianGroup([%d,%d,%d]),SymmetricGroup(2))),chi -> chi[1]))' % (m,m,m)
        degrees=dict(ast.literal_eval(gap.eval(command)))
        q=m**3
        assert degrees=={1:2*q,2:q*(q-1)//2}
        print('GAP',json.dumps({'m':m,'degree_histogram':degrees}),flush=True)
    print('VERSIONS',json.dumps({'python':platform.python_version(),'sage':SAGE_VERSION,'mpmath':mp.__version__,'gap':gap.eval('GAPInfo.Version'),'numberdb':getattr(numberdb,'__version__','unknown')}))
    print('ALL CHECKS PASSED',flush=True)


if __name__=='__main__': main()
