####### 4-16-2007, wdj
### add guava.py to the manual
### fix bug in print and _str_ method - take rref(G) first by default, add NC option
### add zeta_function_fsd_extremal
### formal weight enumerator programs

sage: R0 = PolynomialRing(QQ,"T")      
sage: R = PolynomialRing(CC,"T")      
sage: T = R0.gen()
sage: f = lambda m,i : (1+3*T+2*T^2)^(m)*(1+2*T+2*T^2)^i
sage: c = lambda m,i : f(m,i).list()
sage: q = lambda m,i : [c(m,i)[j]/binomial(4*m+2*i,m+j) for j in range(2*m+1 + 2*i)]
sage: m = 5; Q = R(q(m,2)); rts = Q.roots(); [abs(r) for r in rts]
# zeta poly = R0(q(m,nu))
###### This is for extremal type I codes

sage: f = lambda m,i : (1+3*T+2*T^2)^(m)*(1+2*T+2*T^2)^m*((1+T)^8+14*(1+T)^4*(T)^4+T^8)^i
sage: c = lambda m,i : f(m,i).list()
sage: q = lambda m,i : [c(m,i)[j]/binomial(6*m+8*i,m+j) for j in range(4*m+1 + 8*i)]
sage: m = 5; Q = R(q(m,2)); rts = Q.roots(); [abs(r) for r in rts]
# zeta poly = R0(q(m,nu))/(1-2*T+2*T^2) ??????????????
###### This is for extremal type II codes


sage: f = lambda m,i : (1+3*T+3*T^2)^(m)*((1+T)^4+8*(1+T)*(T)^3)^i
sage: c = lambda m,i : f(m,i).list()
sage: q = lambda m,i : [c(m,i)[j]/binomial(4*m+4*i,m+j) for j in range(2*m+1 + 4*i)]
sage: m = 5; Q = R(q(m,2)); rts = Q.roots(); [abs(r) for r in rts]
# zeta poly = R0(q(m,nu))/(1+3*T^2)  ????????????????
###### This is for extremal type III codes

sage: 
sage: f = lambda m:(1+3*T+2*T^2)^(m)
sage: c = lambda m:f(m).list()
sage: q = lambda m : [c(m)[i]/binomial(4*m,m+i) for i in range(2*m+1)]
sage: m = 5; Q = R(q(m)); rts = Q.roots(); [abs(r) for r in rts]

[0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548]
sage: 
sage: m = 6; Q = R(q(m)); rts = Q.roots(); [abs(r) for r in rts]

[0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548]
sage: m = 7; Q = R(q(m)); rts = Q.roots(); [abs(r) for r in rts]

[0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548]
sage: m = 8; Q = R(q(m)); rts = Q.roots(); [abs(r) for r in rts]

[0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547]
sage: m = 9; Q = R(q(m)); rts = Q.roots(); [abs(r) for r in rts]

[0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186547,
 0.707106781186548,
 0.707106781186547,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186547,
 0.707106781186547,
 0.707106781186547,
 0.707106781186548]
sage: 
sage: f1 = lambda m:(1+3*T+2*T^2)^(m)*(1+2*T+2*T^2)
sage: c1 = lambda m:f1(m).list()
sage: q1 = lambda m : [c1(m)[i]/binomial(4*m+2,m+i) for i in range(2*m+3)]
sage: m = 5; Q = R(q1(m)); rts = Q.roots(); [abs(r) for r in rts]

[0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547]
sage: m = 6; Q = R(q1(m)); rts = Q.roots(); [abs(r) for r in rts]

[0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548]
sage: m = 7; Q = R(q1(m)); rts = Q.roots(); [abs(r) for r in rts]

[0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548]
sage: m = 8; Q = R(q1(m)); rts = Q.roots(); [abs(r) for r in rts]

[0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548]
sage: m = 9; Q = R(q1(m)); rts = Q.roots(); [abs(r) for r in rts]

[0.707106781186547,
 0.707106781186548,
 0.707106781186547,
 0.707106781186547,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186547,
 0.707106781186547,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548]
sage: 
sage: f2 = lambda m:(1+3*T+2*T^2)^(m)*(1+2*T+2*T^2)^2
sage: c2 = lambda m:f2(m).list()
sage: q2 = lambda m : [c2(m)[i]/binomial(4*m+4,m+i) for i in range(2*m+5)]
sage: m = 5; Q = R(q2(m)); rts = Q.roots(); [abs(r) for r in rts]

[0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548]
sage: m = 6; Q = R(q2(m)); rts = Q.roots(); [abs(r) for r in rts]

[0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548]
sage: m = 9; Q = R(q2(m)); rts = Q.roots(); [abs(r) for r in rts]

[0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186547,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548]
sage: 
sage: f4 = lambda m:(1+3*T+2*T^2)^(m)*(1+2*T+2*T^2)^4 # non-existent
sage: c4 = lambda m:f4(m).list()
sage: q4 = lambda m : [c4(m)[i]/binomial(4*m+8,m+i) for i in range(2*m+9)]
sage: m = 5; Q = R(q2(m)); rts = Q.roots(); [abs(r) for r in rts]

[0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548]
sage: m = 15; Q = R(q2(m)); rts = Q.roots(); [abs(r) for r in rts]

[0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186547,
 0.707106781186547,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186547,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548,
 0.707106781186548]

####### formal MDS weight enumerators
##sage: M
##10*y^5 + 30*x*y^4 + 20*x^2*y^3 + 20*x^3*y^2 + x^5
## formal weight enumerator with length 5 and min dist 2

R2 = PolynomialRing(QQ,2,"xy")
x,y = R2.gens()
R2T = PolynomialRing(R2,"T")
T = R2T.gen()
q = 3; m = 5; d = 2
g = sum([T^i for i in range(2*m)])*sum([(q*T)^i for i in range(2*m)])*(y*(1-T)+x*T)^m 
f = g.list()[m-d]
M = (q-1)*f + x^m

q = 2; m = 26; d = 8
g = sum([T^i for i in range(2*m)])*sum([(q*T)^i for i in range(2*m)])*(y*(1-T)+x*T)^m 
f = g.list()[m-d]
M = (q-1)*f + x^m

# to get the formal wt enumerator

coeffs = zeta_poly.list()
Ms = [g.list()[m-d]*(q-1)+x^m for d in range(8,8+13)]
A = sum([coeffs[i]*Ms[i] for i in range(13)])/a0
##sage: A
##y^26 + 650*x^8*y^18 + 845*x^10*y^16 + 2600*x^12*y^14 + 2600*x^14*y^12 + 845*x^16*y^10 + 650*x^18*y^8 + x^26
## for n=26, d=8, q=2

### to get the formal MDS wt enumerator

sage: R = PolynomialRing(QQ,2,"xy")
sage: x,y = R.gens()
sage: f = lambda q,n,m : (x*T+y*(1-T))^(n)*sum([T^i for i in range(m)])*sum([(q*T)^i for i in range(m)])
sage: M = lambda q,n,d,m : (f(q,n,m).list())[d]*(q-1)+x^n