# Independent check of A6-4/A6-6: E6 Bourbaki labelling 1-3-4-5-6 chain, 2 attached to 4.
from fractions import Fraction as Fr
import random
A=[[2,0,-1,0,0,0],[0,2,0,-1,0,0],[-1,0,2,-1,0,0],[0,-1,-1,2,-1,0],[0,0,0,-1,2,-1],[0,0,0,0,-1,2]]
n=6
# (1) weights in fundamental-weight coordinates c_i = <lam, alpha_i^vee>;  s_j(c)_i = c_i - c_j A[j][i]  (A symmetric)
def s_w(j,c): return tuple(c[i]-c[j]*A[j][i] for i in range(n))
def to_dom(c):
    c=tuple(c); word=[]
    while True:
        js=[j for j in range(n) if c[j]<0]
        if not js: return c,word
        j=js[0]; word.append(j+1); c=s_w(j,c)
rho=(1,)*6
lam1=tuple(rho[i]-9*(i==0)-3*(i==5) for i in range(n))   # rho - 9 w1 - 3 w6
lam2=tuple(rho[i]-9*(i==0) for i in range(n))              # rho - 9 w1
d1,w1=to_dom(lam1); d2,w2=to_dom(lam2)
print('lam1*=',lam1,'-> dominant',d1,'word',w1)
print('lam2*=',lam2,'-> dominant',d2,'word',w2)
print('W-conjugate (same dominant rep):',d1==d2,'  rho-omega4 =',tuple(1-(i==3) for i in range(n)))
# (2) rho_L for L = D4 on nodes 2,3,4,5 : positive roots of D4 in simple-root coords (indices 1..4), then pair with alpha_j^vee
def refl_root(j,r): # roots in simple-root coordinates, s_j(r) = r - <r,alpha_j^vee> alpha_j ; <r,alpha_j^vee> = sum_i r_i A[i][j]
    c=sum(r[i]*A[i][j] for i in range(n)); r=list(r); r[j]-=c; return tuple(r)
black=[1,2,3,4]
roots=set(); fr=[tuple(1 if i==j else 0 for i in range(n)) for j in black]
roots|=set(fr)
while fr:
    nf=[]
    for r in fr:
        for j in black:
            r2=refl_root(j,r)
            if r2 not in roots: roots.add(r2); nf.append(r2)
    fr=nf
pos=[r for r in roots if all(x>=0 for x in r)]
assert len(roots)==24 and len(pos)==12
two_rhoL=tuple(sum(r[i] for r in pos) for i in range(n))
pair=lambda r,j: sum(r[i]*A[i][j] for i in range(n))   # <r, alpha_j^vee>
print('2rho_L =',two_rhoL,'  <rho_L,alpha_j^vee> =',[Fr(pair(two_rhoL,j),2) for j in range(n)])
# (3) shifted point: x^{alpha_j^vee} = q^{<rho_L,alpha_j^vee>} * y(phi alpha_j^vee); y=(u,v) on e1=phi a1, e6=phi a6; black j: phi=0
def weight_of(logu,logv):   # lambda* with <lambda*,alpha_j^vee> = log_q x^{alpha_j^vee}
    return tuple(Fr(pair(two_rhoL,j),2)+(logu if j==0 else 0)+(logv if j==5 else 0) for j in range(n))
print('y1=(q^-5,q):  lambda* =',weight_of(-5,1),'  y2=(q^-5,q^4): lambda* =',weight_of(-5,4))
# (4) W_X = S3 on (a,b) coords of Lambda_S : s1(a,b)=(-a,a+b), s6(a,b)=(a+b,-b);  acting on y=(log u, log v) via the dual: (u,v)->(u^-1,uv),(uv,v^-1)
def orb(p):
    s1=lambda p:(-p[0],p[0]+p[1]); s6=lambda p:(p[0]+p[1],-p[1])
    o={p}; fr=[p]
    while fr:
        nf=[]
        for x in fr:
            for y in (s1(x),s6(x)):
                if y not in o: o.add(y); nf.append(y)
        fr=nf
    return o
O1=orb((-5,1)); O2=orb((-5,4))
print('W_X orbit of y1:',sorted(O1)); print('W_X orbit of y2:',sorted(O2)); print('disjoint:',O1.isdisjoint(O2),' sigma(y1)=(1,-5) in orbit of y2:',(1,-5) in O2)
# (5) M_(3,0) at y1,y2:  kappa=(a,b) -> c1=(2a+b)/3, c6=(a+2b)/3 ; y^kappa = u^c1 v^c6
def M(kappa,y):
    vals=[]
    for (a,b) in orb(kappa):
        c1=Fr(2*a+b,3); c6=Fr(a+2*b,3); vals.append(c1*y[0]+c6*y[1])
    return sorted(vals)
for k in [(3,0),(0,3),(1,1),(2,2)]:
    print('M_%s exponents at y1 %s  at y2 %s  equal %s'%(k,M(k,(-5,1)),M(k,(-5,4)),M(k,(-5,1))==M(k,(-5,4))))
# (6) transported W_X action: w1 = w0(D5 on 1..5) w0(D4), w6 = w0(D5 on 2..6) w0(D4) on coweights (coroot coords); check they fix 2rho_L (as coweight via the form) and induce s1,s6 on (a,b)
def refl_cw(j,lam): # coweights in coroot coords: s_j(lam) = lam - <alpha_j,lam> alpha_j^vee, <alpha_j,lam>=sum_i A[j][i] lam_i
    c=sum(A[j][i]*lam[i] for i in range(n)); l=list(lam); l[j]-=c; return tuple(l)
def w0_word(nodes):
    c=[1 if j in nodes else 0 for j in range(n)]; word=[]
    while True:
        js=[j for j in nodes if c[j]>0]
        if not js: return word
        j=js[0]; word.append(j); cj=c[j]; c=[c[k]-cj*A[k][j] for k in range(n)]
def apply(word,lam):
    for j in word: lam=refl_cw(j,lam)
    return lam
wD4=w0_word([1,2,3,4]); wD5a=w0_word([0,1,2,3,4]); wD5b=w0_word([1,2,3,4,5])
W1=lambda lam: apply(wD5a,apply(wD4,lam)); W6=lambda lam: apply(wD5b,apply(wD4,lam))
ab=lambda lam:(sum(A[0][i]*lam[i] for i in range(n)) - 0, sum(A[5][i]*lam[i] for i in range(n)))  # (<alpha_1,phi lam>,<alpha_6,phi lam>) where phi lam = lam - theta lam; use theta
TH=[[-1,0,0,0,0,0],[-1,1,0,0,0,-1],[-2,0,1,0,0,-1],[-2,0,0,1,0,-2],[-1,0,0,0,1,-2],[0,0,0,0,0,-1]]
theta=lambda lam: tuple(sum(TH[i][j]*lam[j] for j in range(n)) for i in range(n))
phi=lambda lam: tuple(lam[i]-theta(lam)[i] for i in range(n))
abphi=lambda lam: ab(phi(lam))
# 2rho_L as a coweight: use the identification roots=coroots (simply laced): coroot coords = simple-root coords
print('w1 fixes 2rho_L:',W1(two_rhoL)==two_rhoL,' w6 fixes 2rho_L:',W6(two_rhoL)==two_rhoL,' lengths',len(wD5a)+len(wD4),len(wD5b)+len(wD4))
ok=True
for _ in range(300):
    lam=tuple(random.randint(-6,6) for _ in range(n)); a,b=abphi(lam)
    if abphi(W1(lam))!=(-a,a+b) or abphi(W6(lam))!=(a+b,-b) or phi(W1(lam))!=W1(phi(lam)) or phi(W6(lam))!=W6(phi(lam)): ok=False
print('w1,w6 induce s1,s6 on (a,b) and commute with phi:',ok)
# (7) phi(alpha_1^vee), phi(alpha_6^vee) in (a,b):
I=[tuple(1 if i==j else 0 for i in range(n)) for j in range(n)]
print('phi(alpha_j^vee) in (a,b):',[abphi(e) for e in I],'  theta^2=1:',all(theta(theta(e))==e for e in I))