import sys; sys.path.insert(0, "/private/tmp/claude-501/-Users-miltonlin-Projects-Trees/75c2c3d9-ef06-4644-b59b-d12583740b0c/scratchpad/solve/A6")
from e6_setup import *
from fractions import Fraction as Fr
from math import gcd
import pickle, random
G=pickle.load(open('e6_points.pkl','rb')); curves=G['curves']
I6=tuple(tuple(1 if i==j else 0 for j in range(n)) for i in range(n))
def apply(w, lam): return tuple(sum(w[j][i]*lam[j] for j in range(n)) for i in range(n))
def wx_images(u,v):
    def s1(p): eu,au,ev,av=p; return (-eu,(-au)%1,eu+ev,(au+av)%1)
    def s6(p): eu,au,ev,av=p; return (eu+ev,(au+av)%1,-ev,(-av)%1)
    orb={(u[0],u[1],v[0],v[1])}; fr=list(orb)
    while fr:
        nf=[]
        for p in fr:
            for q in (s1(p),s6(p)):
                if q not in orb: orb.add(q); nf.append(q)
        fr=nf
    return orb

def inv(w):
    M=[[Fr(w[j][i]) for j in range(n)] for i in range(n)]
    N=[row[:]+[Fr(1 if i==k else 0) for k in range(n)] for i,row in enumerate(M)]
    for c in range(n):
        pp=next(r for r in range(c,n) if N[r][c]!=0); N[c],N[pp]=N[pp],N[c]
        pv=N[c][c]; N[c]=[x/pv for x in N[c]]
        for r in range(n):
            if r!=c and N[r][c]!=0:
                f=N[r][c]; N[r]=[x-f*y for x,y in zip(N[r],N[c])]
    Minv=[[N[i][n+j] for j in range(n)] for i in range(n)]
    return tuple(tuple(int(Minv[i][j]) for i in range(n)) for j in range(n))
random.seed(1)
inconsistent=0; conj=0; nonconj=[]
for w,eqs in curves:
    rows=[(a,b,c) for a,b,c in eqs if (a,b)!=(0,0)]
    if any(c!=0 for a,b,c in eqs if (a,b)==(0,0)): inconsistent+=1; continue
    a0,b0,_=rows[0]; g0=gcd(abs(a0),abs(b0)); p,q=a0//g0,b0//g0
    ks=[]; ms=set()
    for a,b,c in rows:
        k=a//p if p!=0 else b//q
        assert (a,b)==(k*p,k*q)
        ks.append(k); ms.add(Fr(c,k))
    if len(ms)>1: inconsistent+=1; continue
    m=ms.pop(); g=0
    for k in ks: g=gcd(g,abs(k))
    wi=inv(w); cols=[apply(wi,I6[i]) for i in range(n)]
    for j in range(g):           # component: u^p v^q = t^m * e^{2 pi i j/g}
        for trial in range(2):
            # generic point: choose eu, theta_u random; solve for ev, theta_v (q!=0) or vice versa
            if q!=0:
                eu=Fr(random.randint(-50,50),7); ev=(m-p*eu)/q
                au=Fr(random.randint(0,100),101); av=((Fr(j,g)-p*au)/q)%1
            else:
                ev=Fr(random.randint(-50,50),7); eu=(m-q*ev)/p
                av=Fr(random.randint(0,100),101); au=((Fr(j,g)-q*av)/p)%1
            # check all equations
            assert all(a*eu+b*ev==c and ((a*au+b*av)%1)==0 for a,b,c in eqs)
            def xval(c): return (two_rhoL_pair(c)+c[0]*eu+c[5]*ev, (c[0]*au+c[5]*av)%1)
            e1,a1=xval(cols[0]); e6,a6=xval(cols[5]); up=(e1+6,a1); vp=(e6+6,a6)
            if (up[0],up[1],vp[0],vp[1]) in wx_images((eu,au),(ev,av)): conj+=1
            else: nonconj.append((w,(eu,au,ev,av),(up,vp)))
print('rank-1 cases: %d ; inconsistent (no solutions): %d ; generic-point tests W_X-conjugate: %d ; NOT conjugate: %d'%(len(curves),inconsistent,conj,len(nonconj)))