# G7: check the dimension count behind "F_k contains an orbit of dimension >= 7" (verdict page section 2)
# and the fixed-point facts used by the in-run repair of the covering fibre.  Standard library only.
from fractions import Fraction as Fr
from itertools import product
N = 6
edges = [(1,3),(3,4),(4,5),(5,6),(2,4)]
C = [[0]*N for _ in range(N)]
for i in range(N): C[i][i] = 2
for a,b in edges: C[a-1][b-1] = -1; C[b-1][a-1] = -1
def pr(x,y): return sum(x[i]*C[i][j]*y[j] for i in range(N) for j in range(N))
def tofund(x): return [sum(x[i]*C[i][j] for i in range(N)) for j in range(N)]
def inv(M):
    n=len(M); A=[[Fr(M[i][j]) for j in range(n)]+[Fr(int(i==j)) for j in range(n)] for i in range(n)]
    for c in range(n):
        p=next(r for r in range(c,n) if A[r][c]!=0); A[c],A[p]=A[p],A[c]
        pv=A[c][c]; A[c]=[v/pv for v in A[c]]
        for r in range(n):
            if r!=c and A[r][c]!=0:
                f=A[r][c]; A[r]=[A[r][k]-f*A[c][k] for k in range(2*n)]
    return [row[n:] for row in A]
Ci = inv(C)
def fund2root(c): return [sum(Ci[i][j]*c[j] for j in range(N)) for i in range(N)]
# roots
simple=[tuple(int(i==j) for i in range(N)) for j in range(N)]
pos=set(simple); fr=list(simple)
while fr:
    new=[]
    for r in fr:
        f=tofund(list(r))
        for j in range(N):
            if f[j]==-1:
                s=list(r); s[j]+=1; s=tuple(s)
                if s not in pos: pos.add(s); new.append(s)
    fr=new
pos=sorted(pos); roots=pos+[tuple(-v for v in r) for r in pos]
black=[1,2,3,4]
Bv=[[int(i==j) for i in range(N)] for j in black]
G=[[pr(Bv[a],Bv[b]) for b in range(4)] for a in range(4)]; Gi=inv(G)
def Pb(x):
    co=[pr(Bv[a],x) for a in range(4)]
    sol=[sum(Gi[a][b]*co[b] for b in range(4)) for a in range(4)]
    out=[Fr(0)]*N
    for a in range(4):
        for i in range(N): out[i]+=sol[a]*Bv[a][i]
    return out
def theta(x):
    p=Pb(list(x)); return tuple(int(2*p[i]-x[i]) for i in range(N))
def phi(x):
    t=theta(x); return [x[i]-t[i] for i in range(N)]
# roots vanishing on the black torus T' = roots orthogonal to all black coroots
van=[r for r in roots if all(tofund(list(r))[b]==0 for b in black)]
print("roots vanishing on the black torus:", len(van))
e1=phi((1,0,0,0,0,0)); e6=phi((0,0,0,0,0,1))
def kappa(a,b):
    c1=Fr(2*a+b,3); c6=Fr(a+2*b,3)
    return [c1*e1[i]+c6*e6[i] for i in range(N)]
def refl(c,j): return tuple(c[k]-c[j]*C[j][k] for k in range(N))
def orbit(c):
    c=tuple(c); seen={c}; fr=[c]
    while fr:
        new=[]
        for v in fr:
            for j in range(N):
                w=refl(v,j)
                if w not in seen: seen.add(w); new.append(w)
        fr=new
    return seen
for k in (1,3):
    mu=kappa(3*k,0)
    lev=[r for r in roots if pr(list(r),mu)==0]
    nonblack=[r for r in lev if not all(r[i]==0 for i in (0,5))]
    blackr=[r for r in lev if all(r[i]==0 for i in (0,5))]
    pairs=set()
    for r in nonblack:
        t=theta(r); assert t in set(lev)
        pairs.add(frozenset([r,t]))
    selfpaired=[r for r in nonblack if theta(r)==r]
    dimZ=4+len(blackr)+len(pairs)
    print("k=%d: Levi roots of mu_k: %d (black %d, non-black %d in %d theta-pairs, self-paired %d); dim Lie Z_k = %d" %
          (k,len(lev),len(blackr),len(nonblack),len(pairs),len(selfpaired),dimZ))
    lamf=(Fr(3*k-1,2),0,0,0,0,1)
    fixed=[]
    for v in orbit(lamf):
        vr=fund2root(list(v))
        p=phi(vr)
        if all(p[i]==mu[i] for i in range(N)): fixed.append(vr)
    dims=[]
    for nu in fixed:
        nb=sum(1 for r in blackr if pr(list(r),nu)<=0)
        npairs=0
        for P in pairs:
            r,t=tuple(P)
            if pr(list(r),nu)<=0 and pr(list(t),nu)<=0: npairs+=1
        dims.append((4+nb+npairs, dimZ-(4+nb+npairs)))
    print("   fixed points:",len(fixed)," (dim Lie(Z cap P(-nu)), orbit-dim lower bound) per point:", dims)