import sys; sys.path.insert(0,'.')
from lv_module import *
N=12
for n in range(N):
    E=Cs(0,IC[(0,n)]); ICt,mult=kl_reduce(E,(1,n+1)); IC[(1,n+1)]=ICt; IC[(0,n+1)]=Tw(ICt)
# (1) closed form C_(0,n) = v^{n+1} f_n, f_n(O)=max(0,n+1-d(v0,O))
def dv0(O):
    if O in ('P','M'): return 0
    i,m=O; return m if i==0 else m+1
ok=True
for n in range(N):
    pred={}
    for O in ['P','M']+[(i,m) for i in (0,1) for m in range(n+2)]:
        c=max(0,n+1-dv0(O))
        if c: pred[O]={n+1:c}
    if IC[(0,n)]!=pred: ok=False; print("closed form fails",n)
print("closed form C_(0,n)=v^{n+1} f_n:",ok)
# (2) decomposition rule for all n: Cs0 C_(0,n) = C_(1,n+1)+C_(1,n-1); Cs1 C_(0,n)=(v+v^-1)C_(0,n); Cs0 C_(1,n) = (v+1/v) C_(1,n)
ok=True
for n in range(N-1):
    rhs=eadd(IC[(1,n+1)], IC[(1,n-1)] if n>=1 else {})
    ok &= Cs(0,IC[(0,n)])==rhs
    ok &= Cs(1,IC[(0,n)])==escale(IC[(0,n)],{1:1,-1:1})
    ok &= Cs(0,IC[(1,n)])==escale(IC[(1,n)],{1:1,-1:1})
print("decomposition rules hold for n<",N-1,":",ok)
# (3) bar via affine Hecke only: M_{(1,0)} = H_{s0} M_P - v M_M
alt=eadd(Hsbar(0,M('P')), escale(M('M'),mono(-1)),-1)
print("bar(M_B1) via s0 equals via omega:", alt==bar_basis((1,0)))
# (4) Chen-Yi Thm 2(iii) degree property for bar in the 1_O basis: bar(1_v)=q^{-d(v)}1_v + sum b_u 1_u, q^{d(v)} b_u poly in q of degree <= d(v)-d(u)
ok=True
for O in ['P','M']+[(i,m) for i in (0,1) for m in range(8)]:
    b=bar({O:ONE})
    for Op,p in b.items():
        # q = v^-2 ; q^{d(v)} b = v^{-2d(v)} b must be a polynomial in q = v^{-2}: exponents even, <=0, >= -2(d(v)-d(u))
        e=[ex-2*d(O) for ex in p]
        if Op==O:
            if p!={2*d(O):1}: ok=False
            continue
        if any(x%2 or x>0 or x< -2*(d(O)-d(Op)) for x in e): ok=False; print(O,Op,p)
print("Chen-Yi D_delta degree property:",ok)
# (5) function level: annihilator of C_K + C_K' in the dual, on orbits with distance <= L (truncation-stable count)
from fractions import Fraction
L=10
orbs=['P','M']+[(i,m) for i in (0,1) for m in range(L+1)]
idx={O:k for k,O in enumerate(orbs)}
rows=[]
def ind(lst):
    r=[0]*len(orbs)
    for O in lst: r[idx[O]]=1
    return r
rows.append(ind(['P','M',(0,0)])); rows.append(ind(['P','M',(1,0)]))
for n in range(1,L+1):
    rows.append(ind([(0,n),(1,n-1)])); rows.append(ind([(1,n),(0,n-1)]))
# rank over Q
import copy
Mx=[[Fraction(x) for x in r] for r in rows]
rank=0; cols=len(orbs)
for c in range(cols):
    piv=None
    for r in range(rank,len(Mx)):
        if Mx[r][c]!=0: piv=r;break
    if piv is None: continue
    Mx[rank],Mx[piv]=Mx[piv],Mx[rank]
    for r in range(len(Mx)):
        if r!=rank and Mx[r][c]!=0:
            f=Mx[r][c]/Mx[rank][c]; Mx[r]=[a-f*b for a,b in zip(Mx[r],Mx[rank])]
    rank+=1
print("orbits up to distance",L,":",len(orbs)," rank of C_K+C_K' generators:",rank," (the last two orbits (0,L),(1,L) are only partially covered)")
