# Version for filtered Hall bases.
# The only difference is the use of files bas_comm.Nf 
# instead of bas_comm.N

# Expand the log of the KZ Drinfeld associator up to a certain degree
# over the Hall bases
# Needs: mzv-1-12.txt (the French tables)
#        bas_comm.2f -- bas_comm.12f
#        assolist.2 -- 12
# no subtitution for zeta(2)!

interface(prettyprint=0);
interface(quiet=true);
interface(echo=0);
interface(screenwidth=255);
printlevel:=0;
#writeto(`asso_log12f.res`);

with(linalg):
read(`mzv3-12.mpl`):
read(`assolist.2-8`):
read(`bas_comm.2-8f`):

# Concatenation multiplication:
cm:= proc(p1,p2)
local e1,e2,t1,t2,v1,v2,n1,n2,i1,i2,c1,c2,v,r;
e1:=expand(p1);
e2:=expand(p2);
r:=0;
for i1 from 1 to nops(e1) do
for i2 from 1 to nops(e2) do
t1:=op(i1,e1);
t2:=op(i2,e2);
n1:=nops(t1);
n2:=nops(t2);
v1:=op(n1,t1);
v2:=op(n2,t2);
c1:=coeff(t1,v1);
c2:=coeff(t2,v2);
r:=r+c1*c2*cat(v1,v2);
od;od;
RETURN(r);
end:

# assign A->0, B->1
digit01:=proc(a)
if a=65 then 0 else 1 fi;
end:

# assign a number to a monomial in letters A,B
numvar:= proc(a)
local la,i,r;
la:=convert(a,bytes);
r:=0;
for i from 1 to nops(la) do
r:=r*2+digit01(la[i]); od;
RETURN(r);
end:

# Convert a linear combination of variables
# AB,BA,AAB,ABA,ABB,BAA,BAB,BBA,AAAB... into a vector
tovect:=proc(n,p)
local a,k,i,c,L,v;
L:=2^n-2;
for i from 1 to L do v[i]:=0 od;
for i from 1 to nops(p) do
a:=op(i,p);
if nops(a)=1 then c:=1;k:=numvar(a)
else c:=op(1,a);k:=numvar(op(2,a)) fi;
v[k]:=v[k]+c;
od;
RETURN([seq(v[i],i=1..L)]);
end:

# extract the coefficient of v in p:
coeffm:= proc(p,v)
local t,v1,v2,n,r;
t:=whattype(v);
if t='symbol' then r:=coeff(p,v)
elif t='function' then r:=coeff(p,v)
elif t='`*`' then 
v1:=op(1,v); v2:=product(op(ii,v),ii=2..nops(v)); r:=coeffm(coeff(p,v1),v2)
elif t='`^`' then 
v1:=op(1,v); n:=op(2,v); r:=coeff(p,v1,n)
fi;
RETURN(r);
end:

# degree 2

log2:=ass2:
matr2t:=transpose(matr2):
listz:=[z2]:
for i from 1 to 1 do
log2z||i:=simplify(coeffm(log2,listz[i])):
end:
# uncomment for debugging:
# simplify(log2-sum('log2z||j'*listz[j],j=1..1));

for i from 1 to 1 do
n:=numer(log2z||i):
d:=denom(log2z||i):
coef[i]:=listz[i]/d:
nv:=tovect(2,n):
v[i]:=linsolve(matr2t,nv):
od:

len:=1:
lv:=1:
for i from 1 to lv do
print(simplify(sum(v[j][i]*coef[j],j=1..len)));
od;

# degree 3
log3:=ass3:
matr3t:=transpose(matr3):
listz:=[z3]:
for i from 1 to 1 do
log3z||i:=simplify(coeffm(log3,listz[i])):
end:
# uncomment for debugging:
# simplify(log3-sum('log3z||j'*listz[j],j=1..1));

len:=1:
for i from 1 to 1 do
n:=numer(log3z||i):
d:=denom(log3z||i):
coef[i]:=listz[i]/d:
nv:=tovect(3,n):
v[i]:=linsolve(matr3t,nv):
od:

lv:=2:
for i from 1 to lv do
print(simplify(sum(v[j][i]*coef[j],j=1..len)));
od;

# degree 4
log4:=ass4-1/2*cm(ass2,ass2):
matr4t:=transpose(matr4):
listz:=[z2^2]:
for i from 1 to 1 do
log4z||i:=simplify(coeffm(log4,listz[i])):
end:
# uncomment for debugging:
# simplify(log4-sum('log4z||j'*listz[j],j=1..1));

for i from 1 to 1 do
n:=numer(log4z||i):
d:=denom(log4z||i):
coef[i]:=listz[i]/d:
nv:=tovect(4,n):
v[i]:=linsolve(matr4t,nv):
od:

lv:=3:
for i from 1 to lv do
print(simplify(sum(v[j][i]*coef[j],j=1..len)));
od;


# degree 5
log5:=ass5-1/2*(cm(ass2,ass3)+cm(ass3,ass2)):
matr5t:=transpose(matr5):
listz:=[z2*z3,z5]:
len:=nops(listz):
for i from 1 to len do
log5z||i:=simplify(coeffm(log5,listz[i])):
end:
# uncomment for debugging:
# simplify(log5-sum('log5z||j'*listz[j],j=1..len));

for i from 1 to len do
n:=numer(log5z||i):
d:=denom(log5z||i):
coef[i]:=listz[i]/d:
nv:=tovect(5,n):
v[i]:=linsolve(matr5t,nv):
od:

# now, using the collected data, write the answer expanded over the
# basis commutators, not over the MZV's:

lv:=6:
for i from 1 to lv do
print(simplify(sum(v[j][i]*coef[j],j=1..len)));
od;

# degree 6
log6:=ass6-1/2*(cm(ass2,ass4)+cm(ass3,ass3)+cm(ass4,ass2))
      +1/3*cm(cm(ass2,ass2),ass2):
matr6t:=transpose(matr6):
listz:=[z2^3, z3^2]:
len:=nops(listz):
for i from 1 to len do
log6z||i:=simplify(coeffm(log6,listz[i])):
end:
# uncomment for debugging:
# simplify(log6-sum('log6z||j'*listz[j],j=1..len));

for i from 1 to len do
n:=numer(log6z||i):
d:=denom(log6z||i):
coef[i]:=listz[i]/d:
nv:=tovect(6,n):
v[i]:=linsolve(matr6t,nv):
od:

lv:=9:
for i from 1 to lv do
print(simplify(sum(v[j][i]*coef[j],j=1..len)));
od;

# degree 7
log7:=ass7-1/2*(cm(ass2,ass5)+cm(ass3,ass4)+cm(ass4,ass3)+cm(ass5,ass2))
+1/3*(cm(cm(ass2,ass2),ass3)+cm(ass2,cm(ass3,ass2))+cm(ass3,cm(ass2,ass2))):
listz:=[z2^2*z3,z2*z5,z7]:
matr7t:=transpose(matr7):
len:=nops(listz):
for i from 1 to len do
log7z||i:=simplify(coeffm(log7,listz[i])):
end:
# uncomment for debugging:
# simplify(log7-sum('log7z||j'*listz[j],j=1..len));

for i from 1 to len do
n:=numer(log7z||i):
d:=denom(log7z||i):
coef[i]:=listz[i]/d:
nv:=tovect(7,n):
v[i]:=linsolve(matr7t,nv):
od:

lv:=18:
for i from 1 to lv do
print(simplify(sum(v[j][i]*coef[j],j=1..len)));
od;

# degree 8

log8:=ass8-1/2*(cm(ass2,ass6)+cm(ass3,ass5)+cm(ass4,ass4)+cm(ass5,ass3)
+cm(ass6,ass2))+1/3*(cm(cm(ass2,ass2),ass4)+cm(ass2,cm(ass4,ass2))
+cm(ass4,cm(ass2,ass2))+cm(cm(ass2,ass3),ass3)+cm(ass3,cm(ass2,ass3))
+cm(ass3,cm(ass3,ass2)))-1/4*cm(ass2,cm(ass2,cm(ass2,ass2))):
listz:=[z2^4, z2*z3^2, z3*z5, z62]:
matr8t:=transpose(matr8):
len:=nops(listz):
for i from 1 to len do
log8z||i:=simplify(coeffm(log8,listz[i])):
end:
# uncomment for debugging:
# simplify(log8-sum('log8z||j'*listz[j],j=1..len));

for i from 1 to len do
n:=numer(log8z||i):
d:=denom(log8z||i):
coef[i]:=listz[i]/d:
nv:=tovect(8,n):
v[i]:=linsolve(matr8t,nv):
od:

lv:=30:
for i from 1 to lv do
print(simplify(sum(v[j][i]*coef[j],j=1..len)));
od;

# degree 9
read(`assolist.09`):
log9:=ass9-1/2*(cm(ass2,ass7)+cm(ass3,ass6)+cm(ass4,ass5)+cm(ass5,ass4)
+cm(ass6,ass3)+cm(ass7,ass2))
+1/3*(cm(cm(ass2,ass2),ass5)+cm(ass2,cm(ass5,ass2))
+cm(ass5,cm(ass2,ass2))+cm(ass3,cm(ass3,ass3))
+cm(cm(ass2,ass3),ass4)+cm(cm(ass2,ass4),ass3)+cm(cm(ass3,ass2),ass4)
+cm(cm(ass3,ass4),ass2)+cm(cm(ass4,ass2),ass3)+cm(cm(ass4,ass3),ass2))
-1/4*(cm(ass2,cm(ass2,cm(ass2,ass3)))+cm(ass2,cm(ass2,cm(ass3,ass2)))
+cm(ass2,cm(ass3,cm(ass2,ass2)))+cm(ass3,cm(ass2,cm(ass2,ass2)))):

listz:=[z2^3*z3, z3^3, z2^2*z5, z2*z7, z9]:

read(`bas_comm.9f`):
matr9t:=transpose(matr9):
len:=nops(listz):
for i from 1 to len do
log9z||i:=simplify(coeffm(log9,listz[i])):
end:
# comment/uncomment for debugging:
difference:=simplify(log9-sum('log9z||j'*listz[j],j=1..len)):
if difference=0 then print('difference=0') fi;

for i from 1 to len do
n:=numer(log9z||i):
d:=denom(log9z||i):
coef[i]:=listz[i]/d:
nv:=tovect(9,n):
v[i]:=linsolve(matr9t,nv):
od:

lv:=56:
for i from 1 to lv do
print(simplify(sum(v[j][i]*coef[j],j=1..len)));
od;

# degree 10

read(`assolist.10`):
logphi10:=ass10-1/2*(cm(ass2,ass8)+cm(ass3,ass7)+cm(ass4,ass6)+cm(ass5,ass5)
+cm(ass6,ass4)+cm(ass7,ass3)+cm(ass8,ass2))
+1/3*(cm(cm(ass2,ass2),ass6)+cm(ass2,cm(ass6,ass2))+cm(ass6,cm(ass2,ass2))
+cm(ass3,cm(ass3,ass4))+cm(ass3,cm(ass4,ass3))+cm(ass4,cm(ass3,ass3))
+cm(ass2,cm(ass4,ass4))+cm(ass4,cm(ass2,ass4))+cm(ass4,cm(ass4,ass2))
+cm(cm(ass2,ass3),ass5)+cm(cm(ass2,ass5),ass3)+cm(cm(ass3,ass2),ass5)
+cm(cm(ass3,ass5),ass2)+cm(cm(ass5,ass2),ass3)+cm(cm(ass5,ass3),ass2))
-1/4*(cm(ass2,cm(ass2,cm(ass2,ass4)))+cm(ass2,cm(ass2,cm(ass4,ass2)))
+cm(ass2,cm(ass4,cm(ass2,ass2)))+cm(ass4,cm(ass2,cm(ass2,ass2)))
+cm(ass2,cm(ass2,cm(ass3,ass3)))+cm(ass2,cm(ass3,cm(ass2,ass3))) 
+cm(ass2,cm(ass3,cm(ass3,ass2)))+cm(ass3,cm(ass2,cm(ass2,ass3)))
+cm(ass3,cm(ass2,cm(ass3,ass2)))+cm(ass3,cm(ass3,cm(ass2,ass2))))
+1/5*cm(ass2,cm(ass2,cm(ass2,cm(ass2,ass2)))):

listz:=[z2^5, z2^2*z3^2, z2*z3*z5, z2*z62, z3*z7, z5^2, z82]:
len:=nops(listz):

read(`bas_comm.10f`):
matr10t:=transpose(matr10):

for i from 1 to len do
log10z||i:=simplify(coeffm(logphi10,listz[i])):
od:

# comment/uncomment for debugging:
difference:=simplify(logphi10-sum('log10z||j'*listz[j],j=1..len)):
if difference=0 then print('difference=0') fi;

for i from 1 to len do
n:=numer(log10z||i):
d:=denom(log10z||i):
coef[i]:=listz[i]/d:
nv:=tovect(10,n):
v[i]:=linsolve(matr10t,nv):
od:

lv:=99:
for i from 1 to lv do
print(simplify(sum(v[j][i]*coef[j],j=1..len)));
od;

# degree 11

read(`assolist.11`):

logphi11:=ass11

-1/2*(cm(ass2,ass9)+cm(ass3,ass8)+cm(ass4,ass7)+cm(ass5,ass6)
+cm(ass6,ass5)+cm(ass7,ass4)+cm(ass8,ass3)+cm(ass9,ass2))

+1/3*(
cm(cm(ass2,ass2),ass7)+cm(ass2,cm(ass7,ass2))+cm(ass7,cm(ass2,ass2))

+cm(cm(ass2,ass3),ass6)+cm(cm(ass2,ass6),ass3)+cm(cm(ass3,ass2),ass6)
+cm(cm(ass3,ass6),ass2)+cm(cm(ass6,ass2),ass3)+cm(cm(ass6,ass3),ass2)

+cm(cm(ass2,ass4),ass5)+cm(cm(ass2,ass5),ass4)+cm(cm(ass4,ass2),ass5)
+cm(cm(ass4,ass5),ass2)+cm(cm(ass5,ass2),ass4)+cm(cm(ass5,ass4),ass2)

+cm(ass3,cm(ass3,ass5))+cm(ass3,cm(ass5,ass3))+cm(ass5,cm(ass3,ass3))

+cm(ass3,cm(ass4,ass4))+cm(ass4,cm(ass3,ass4))+cm(ass4,cm(ass4,ass3))
)

-1/4*(
cm(ass2,cm(ass2,cm(ass2,ass5)))+cm(ass2,cm(ass2,cm(ass5,ass2)))
+cm(ass2,cm(ass5,cm(ass2,ass2)))+cm(ass5,cm(ass2,cm(ass2,ass2)))

+cm(ass2,cm(ass2,cm(ass3,ass4)))+cm(ass2,cm(ass3,cm(ass2,ass4))) 
+cm(ass2,cm(ass3,cm(ass4,ass2)))+cm(ass3,cm(ass2,cm(ass2,ass4)))
+cm(ass3,cm(ass2,cm(ass4,ass2)))+cm(ass3,cm(ass4,cm(ass2,ass2)))
+cm(ass2,cm(ass2,cm(ass4,ass3)))+cm(ass2,cm(ass4,cm(ass2,ass3))) 
+cm(ass2,cm(ass4,cm(ass3,ass2)))+cm(ass4,cm(ass2,cm(ass2,ass3)))
+cm(ass4,cm(ass2,cm(ass3,ass2)))+cm(ass4,cm(ass3,cm(ass2,ass2)))

+cm(ass3,cm(ass3,cm(ass3,ass2)))+cm(ass3,cm(ass3,cm(ass2,ass3)))
+cm(ass3,cm(ass2,cm(ass3,ass3)))+cm(ass2,cm(ass3,cm(ass3,ass3)))
)

+1/5*(
cm(ass2,cm(ass2,cm(ass2,cm(ass2,ass3))))
+cm(ass2,cm(ass2,cm(ass2,cm(ass3,ass2))))
+cm(ass2,cm(ass2,cm(ass3,cm(ass2,ass2))))
+cm(ass2,cm(ass3,cm(ass2,cm(ass2,ass2))))
+cm(ass3,cm(ass2,cm(ass2,cm(ass2,ass2))))
)
:

listz:=[z2^4*z3, z2*z3^3, z2^3*z5, z3^2*z5, 
      z3*z62, z2^2*z7, z821, z2*z9, zB]:
len:=nops(listz):


for i from 1 to len do
log11z||i:=simplify(coeffm(logphi11,listz[i])):
od:

difference:=simplify(logphi11-sum('log11z||j'*listz[j],j=1..len)):
if difference=0 then print('difference=0') fi;

read(`bas_comm.11f`):
matr11t:=transpose(matr11):

for i from 1 to len do
n:=numer(log11z||i):
d:=denom(log11z||i):
coef[i]:=listz[i]/d:
nv:=tovect(11,n):
v[i]:=linsolve(matr11t,nv):
od:

lv:=186:
for i from 1 to lv do
print(simplify(sum(v[j][i]*coef[j],j=1..len)));
od;


# degree 12

read(`assolist.12`):

logphi12:=ass12

-1/2*(cm(ass2,ass10)+cm(ass3,ass9)+cm(ass4,ass8)+cm(ass5,ass7)
+cm(ass6,ass6)+cm(ass7,ass5)+cm(ass8,ass4)+cm(ass9,ass3)+cm(ass10,ass2))

+1/3*(
cm(cm(ass2,ass2),ass8)+cm(ass2,cm(ass8,ass2))+cm(ass8,cm(ass2,ass2))

+cm(cm(ass2,ass3),ass7)+cm(cm(ass2,ass7),ass3)+cm(cm(ass3,ass2),ass7)
+cm(cm(ass3,ass7),ass2)+cm(cm(ass7,ass2),ass3)+cm(cm(ass7,ass3),ass2)

+cm(cm(ass2,ass4),ass6)+cm(cm(ass2,ass6),ass4)+cm(cm(ass4,ass2),ass6)
+cm(cm(ass4,ass6),ass2)+cm(cm(ass6,ass2),ass4)+cm(cm(ass6,ass4),ass2)

+cm(ass2,cm(ass5,ass5))+cm(ass5,cm(ass2,ass5))+cm(ass5,cm(ass5,ass2))

+cm(ass3,cm(ass3,ass6))+cm(ass3,cm(ass6,ass3))+cm(ass6,cm(ass3,ass3))

+cm(ass3,cm(ass4,ass5))+cm(ass3,cm(ass5,ass4))+cm(ass4,cm(ass3,ass5))
+cm(ass4,cm(ass5,ass3))+cm(ass5,cm(ass3,ass4))+cm(ass5,cm(ass4,ass3))

+cm(ass4,cm(ass4,ass4))
)

-1/4*(
cm(ass2,cm(ass2,cm(ass2,ass6)))+cm(ass2,cm(ass2,cm(ass6,ass2)))
+cm(ass2,cm(ass6,cm(ass2,ass2)))+cm(ass6,cm(ass2,cm(ass2,ass2)))

+cm(ass2,cm(ass2,cm(ass3,ass5)))+cm(ass2,cm(ass3,cm(ass2,ass5))) 
+cm(ass2,cm(ass3,cm(ass5,ass2)))+cm(ass3,cm(ass2,cm(ass2,ass5)))
+cm(ass3,cm(ass2,cm(ass5,ass2)))+cm(ass3,cm(ass5,cm(ass2,ass2)))
+cm(ass2,cm(ass2,cm(ass5,ass3)))+cm(ass2,cm(ass5,cm(ass2,ass3))) 
+cm(ass2,cm(ass5,cm(ass3,ass2)))+cm(ass5,cm(ass2,cm(ass2,ass3)))
+cm(ass5,cm(ass2,cm(ass3,ass2)))+cm(ass5,cm(ass3,cm(ass2,ass2)))

+cm(ass2,cm(ass2,cm(ass4,ass4)))+cm(ass2,cm(ass4,cm(ass2,ass4)))
+cm(ass2,cm(ass4,cm(ass4,ass2)))+cm(ass4,cm(ass2,cm(ass2,ass4)))
+cm(ass4,cm(ass2,cm(ass4,ass2)))+cm(ass4,cm(ass4,cm(ass2,ass2)))

+cm(ass2,cm(ass3,cm(ass3,ass4)))+cm(ass2,cm(ass3,cm(ass4,ass3)))
+cm(ass2,cm(ass4,cm(ass3,ass3)))+cm(ass3,cm(ass2,cm(ass3,ass4)))
+cm(ass3,cm(ass2,cm(ass4,ass3)))+cm(ass4,cm(ass2,cm(ass3,ass3)))
+cm(ass3,cm(ass3,cm(ass2,ass4)))+cm(ass3,cm(ass4,cm(ass2,ass3)))
+cm(ass4,cm(ass3,cm(ass2,ass3)))+cm(ass3,cm(ass3,cm(ass4,ass2)))
+cm(ass3,cm(ass4,cm(ass3,ass2)))+cm(ass4,cm(ass3,cm(ass3,ass2)))

+cm(ass3,cm(ass3,cm(ass3,ass3)))
)

+1/5*(
cm(ass2,cm(ass2,cm(ass2,cm(ass2,ass4))))
+cm(ass2,cm(ass2,cm(ass2,cm(ass4,ass2))))
+cm(ass2,cm(ass2,cm(ass4,cm(ass2,ass2))))
+cm(ass2,cm(ass4,cm(ass2,cm(ass2,ass2))))
+cm(ass4,cm(ass2,cm(ass2,cm(ass2,ass2))))

+cm(ass2,cm(ass2,cm(ass2,cm(ass3,ass3))))
+cm(ass2,cm(ass2,cm(ass3,cm(ass2,ass3))))
+cm(ass2,cm(ass3,cm(ass2,cm(ass2,ass3))))
+cm(ass3,cm(ass2,cm(ass2,cm(ass2,ass3))))
+cm(ass2,cm(ass2,cm(ass3,cm(ass3,ass2))))
+cm(ass2,cm(ass3,cm(ass2,cm(ass3,ass2))))
+cm(ass3,cm(ass2,cm(ass2,cm(ass3,ass2))))
+cm(ass2,cm(ass3,cm(ass3,cm(ass2,ass2))))
+cm(ass3,cm(ass2,cm(ass3,cm(ass2,ass2))))
+cm(ass3,cm(ass3,cm(ass2,cm(ass2,ass2))))
)

-1/6*cm(ass2,cm(ass2,cm(ass2,cm(ass2,cm(ass2,ass2)))))
:
listz:=[z2^6, z2^3*z3^2, z3^4, z2^2*z3*z5, z2*z5^2, z2^2*z62, 
      z2*z3*z7, z5*z7, z3*z9, z2*z82, z8211, zA2]:
len:=nops(listz):

for i from 1 to len do
log12z||i:=simplify(coeffm(logphi12,listz[i])):
od:

difference:=simplify(logphi12-sum('log12z||j'*listz[j],j=1..len)):
if difference=0 then print('difference=0') fi;

read(`bas_comm.12f`):
matr12t:=transpose(matr12):

for i from 1 to len do
n:=numer(log12z||i):
d:=denom(log12z||i):
coef[i]:=listz[i]/d:
nv:=tovect(12,n):
v[i]:=linsolve(matr12t,nv):
od:

lv:=335:
for i from 1 to lv do
print(simplify(sum(v[j][i]*coef[j],j=1..len)));
od;

quit;