#OK to post homework
#Dayoon Kim, 2024-2-4, HW5

#Old code C4, C5:
Help4:=proc(): print(`Fqn(q,n), HD(u,v), RV(q,n) , RC(q,n,d,K), SPB(q,n,t)`):end:
Help5:=proc(): print(`Nei(q,c), SP(q,c,t), GRC(q,n,d), MinD(C), CV(S,n), BDtoC(BD,n)`):end:
Help:=proc(): print(`Info(q, d, n_min, n_max), PlotkinBound(n,d)`):end:

#Alphabet {0,1,...q-1}, Fqn(q,n): {0,1,...,q-1}^n
Fqn:=proc(q,n) local S,a,v:
option remember:
if n=0 then
 RETURN({[]}):
fi:
S:=Fqn(q,n-1):
{seq(seq([op(v),a],a=0..q-1), v in S)}:
end:

#HD(u,v): The Hamming distance between two words (of the same length)
HD:=proc(u,v) local i,co:
co:=0:
for i from 1 to nops(u) do
  if u[i]<>v[i] then
      co:=co+1:
  fi:
od:
co:
end:

#SPB(q,n,d): The best you can hope for (as far as the size of C) for q-ary (n,2*t+1) code
SPB:=proc(q,n,t) local i:
trunc(q^n/add(binomial(n,i)*(q-1)^i,i=0..t)):
end:


#RV(q,n): A random word of length n in {0,1,..,q-1}
RV:=proc(q,n) local ra,i:
ra:=rand(0..q-1):
[seq( ra(), i=1..n)]:
end:

#RC(q,n,d,K): inputs q,n,d, and K and keeps picking K+1 (thanks to Nuray) random vectors
#whenver the new vector is not distance <=d-1 from all the previous one 
RC:=proc(q,n,d,K) local C,v,i,c:
C:={RV(q,n)}:
for i from 1 to K do
 v:=RV(q,n):
  if min(seq(HD(v,c), c in C))>=d then
   C:=C union {v}:
  fi:
od:
C:
end:

#Nei(q,c): all the neighbors of the vector c in Fqn(q,n)
Nei:=proc(q,c) local n,i,a: 
n:=nops(c):
{seq(seq([op(1..i-1,c),a,op(i+1..n,c)], a=0..q-1) , i=1..n)}:
end:

#SP(q,c,t): the set of all vectors in Fqn(q,n) whose distance is <=t from c
SP:=proc(q,c,t) local S,s,i:
S:={c}:
for i from 1 to t do
 S:=S union {seq(op(Nei(q,s)),s in S)}:
od:
S:
end:

#GRC(q,n,d):
GRC:=proc(q,n,d) local S,A,v:
A:=Fqn(q,n):
S:={}:
while A<>{} do:
 v:=A[1]:
 S:=S union {v}:
 A:=A minus SP(q,v,d-1):
od:
S:
end:

#MinD(C): The minimal (Hamming) distance of the code C
MinD:=proc(C) local i,j:
min( seq(seq(HD(C[i],C[j]),j=i+1..nops(C)), i=1..nops(C))):
end:

#CV(S,n): the characteristic vector of the subset S of {1,...,n}
CV:=proc(S,n) local v,i:
v:=[]:
for i from 1 to n do
 if member(i,S) then
  v:=[op(v),1]:
 else
  v:=[op(v),0]:
 fi:
od:
v:
end:

BDtoC:=proc(BD,n) local s, C:
C:={seq(CV(s,n),s in BD)}:
C:=C union subs({0=1,1=0},C):
C union {[0$n],[1$n]}:
end:

#End of Old code C5

#-------------------------------------------
#HW5-3.
#HW5-4.
#HW5-5.
Info := proc(q, d, n_min, n_max) local n,C,M,SPB_val:
    for n from n_min to n_max do
        C := GRC(q, n, d);
        M := nops(C);
        SPB_val := SPB(q, n, (d-1)/2);
        
        print("q =", q, ", n =", n, ", d =", d);
        print("Cardinality of the code (M):", M);
        print("Sphere Packing Bound (SPB):", SPB_val);
        print("Difference (SPB - M):", SPB_val - M);
        print("");
    od:
end:

#-------------------------------------------
#HW5-7. PlotkinBound(n,d)

# Define the Plotkin Bound function
PlotkinBound:=proc(n,d) local q, M:
  q := 2;
  if n <= 2*d then
    M := 2*trunc(d/abs(d-n));
  else
    M := n;
  fi;
  M;
end:

# Compare the Sphere Packing Bound and the Plotkin Bound
CompareSPBnPB:=proc() local d,n,SPB_value,Plotkin_value:
for d from 2 to 20 do
  for n from 2 to 2*d do
    SPB_value := SPB(2, n, d);
    Plotkin_value := PlotkinBound(n, d);
    if SPB_value > Plotkin_value then
      printf("For n=%d and d=%d, the Sphere Packing Bound (%d) is larger than the Plotkin Bound (%d).\n", n, d, SPB_value, Plotkin_value);
    elif SPB_value < Plotkin_value then
      printf("For n=%d and d=%d, the Plotkin Bound (%d) is larger than the Sphere Packing Bound (%d).\n", n, d, Plotkin_value, SPB_value);
    else
      printf("For n=%d and d=%d, the Sphere Packing Bound and the Plotkin Bound are equal (%d).\n", n, d, SPB_value);
    fi;
  od;
od;