"""Lemma 23.3, step 2: (U_{Wh,c}, lambda_c) ~ (U_alt, its character), at the level of Lie algebras in so_{2m+1},
m=2n, for several rational values of c: row separation (Lemma 22.3 with r=n, b=c/2), interleaving
(Lemma 22.4 with A_i=e_i, B_j=e_{n+j}, d=1), then the isometry A_i->e_{2i-1}, B_j->e_{2j}, e_m<->f_m, u->-u.
Every hypothesis of the exchange lemma is checked on the Lie algebra (closure, invariance of the character, pairing)."""
import sympy as sp, itertools
def run(n,c):
    m=2*n; N=2*m+1
    E_=lambda i:i-1; F_=lambda i:m+i-1; U=2*m
    def Em(a,b):
        M=sp.zeros(N,N); M[a,b]=1; return M
    def rv(x):
        idx=[(i+1,t) for i,t in enumerate(x) if t]
        if len(idx)==2:
            (i,a),(j,b)=idx
            if a==1 and b==-1: return Em(E_(i),E_(j))-Em(F_(j),F_(i))
            if a==1 and b==1: return Em(E_(i),F_(j))-Em(E_(j),F_(i))
            if a==-1 and b==1: return Em(E_(j),E_(i))-Em(F_(i),F_(j))
        (i,a),=idx
        return 2*Em(E_(i),U)-Em(U,F_(i))
    def v(*t):
        x=[0]*m
        for sg,i in t: x[i-1]+=sg
        return tuple(x)
    b=sp.Rational(c)/2
    L=lambda X,i: X[E_(i),E_(m)]+b*X[E_(i),F_(m)]
    def lam_c(X):
        return sum(X[E_(i),E_(i+1)] for i in range(1,n))+L(X,n)+sum(X[E_(i),E_(i+1)] for i in range(n+1,m-1))+X[E_(m-1),F_(m)]
    # basis of u_{Wh,c}
    basis=[]
    for i in range(1,m+1):
        if i!=m: basis.append(rv(v((1,i))))
        for j in range(i+1,m+1):
            if j==m and n<i<m: continue
            basis.append(rv(v((1,i),(-1,j)))); basis.append(rv(v((1,i),(1,j))))
    for i in range(n+1,m):
        basis.append(-b*rv(v((1,i),(-1,m)))+rv(v((1,i),(1,m))))
    def vecs(B): return sp.Matrix([list(x) for x in B]).T
    def inspan(B,X):
        A=vecs(B); return A.rank()==A.row_join(sp.Matrix(list(X))).rank()
    br=lambda X,Y:X*Y-Y*X
    log=[]
    def chk_alg(B,zeta,what):
        A=vecs(B); rk=A.rank()
        for X,Y in itertools.combinations(B,2):
            Z=br(X,Y)
            if Z!=sp.zeros(N,N):
                if A.row_join(sp.Matrix(list(Z))).rank()!=rk: log.append(what+': not closed'); return
                if zeta(Z)!=0: log.append(what+': character on commutator'); return
    def exchange(B,zeta,x,y,what):
        Uc=[t for t in B if t!=x]
        assert len(Uc)==len(B)-1,what
        chk_alg(Uc,zeta,what+' common')
        if zeta(x)!=0 or zeta(y)!=0: log.append(what+': character on X or Y')
        Z=br(x,y)
        if not inspan(Uc,Z) or zeta(Z)==0: log.append(what+': pairing')
        A=vecs(Uc); rk=A.rank()
        for nm,w in (('X',x),('Y',y)):
            for t in Uc:
                Z=br(w,t)
                if Z!=sp.zeros(N,N):
                    if A.row_join(sp.Matrix(list(Z))).rank()!=rk: log.append(f'{what}: {nm} does not normalize')
                    elif zeta(Z)!=0: log.append(f'{what}: {nm} changes character')
        return Uc+[y]
    chk_alg(basis,lam_c,'initial')
    for j in range(n+1,m):
        basis=exchange(basis,lam_c,rv(v((1,n),(-1,j))),rv(v((1,j),(-1,m))),f'rowsep j={j}')
    A_=lambda i:i; B_=lambda j:n+j; d=1
    for i in range(n-1,0,-1):
        for j in range(1,i-d+1):
            basis=exchange(basis,lam_c,rv(v((1,A_(i)),(-1,B_(j)))),rv(v((1,B_(j)),(-1,A_(i+1)))),f'interleave ({i},{j})')
    # the isometry
    P=sp.zeros(N,N)
    for i in range(1,n+1): P[E_(2*i-1),E_(i)]=1; P[F_(2*i-1),F_(i)]=1
    for j in range(1,n): P[E_(2*j),E_(n+j)]=1; P[F_(2*j),F_(n+j)]=1
    P[F_(m),E_(m)]=1; P[E_(m),F_(m)]=1; P[U,U]=-1
    Q=sp.zeros(N,N)
    for i in range(1,m+1): Q[E_(i),F_(i)]=Q[F_(i),E_(i)]=1
    Q[U,U]=2
    iso=(P.T*Q*P==Q) and P.det()==1
    # U_alt
    alt=[]
    for i in range(1,m+1):
        if i!=m: alt.append(rv(v((1,i))))
        for j in range(i+1,m+1):
            if not (j==i+1 and i%2==0): alt.append(rv(v((1,i),(-1,j))))
            alt.append(rv(v((1,i),(1,j))))
    def zeta_alt(X):
        return sum(X[E_(i),E_(i+2)] for i in range(1,m-1))+X[E_(m-1),F_(m)]+b*X[E_(m-1),E_(m)]
    Pi=P.inv()
    moved=[P*X*Pi for X in basis]
    same_space=len(moved)==len(alt) and vecs(moved).rank()==vecs(alt).rank()==vecs(moved+alt).rank()
    same_char=all(sp.simplify(zeta_alt(P*X*Pi)-lam_c(X))==0 for X in basis)
    chk_alg(alt,zeta_alt,'U_alt')
    print(f'n={n} c={c}: problems={len(log)} isometry det 1={iso} image = U_alt: {same_space} characters agree: {same_char}',log[:3])
import sys
for n in (2,3,4):
    for c in ((0,2,sp.Rational(-4,3)) if n<4 else (2,)): run(n,c)