######################################################################
## Marko.txt Save this file as   Marko.txt to use it,                #
# stay in the                                                        #
## same directory, get into Maple (by typing: maple <Enter> )        #
## and then type:  read `Marko.txt` <Enter>                          #
## Then follow the instructions given there                          #
##                                                                   #
## Written by Doron Zeilberger, Rutgers University ,                 #
##  DoronZeil at gmail dot com                                       # 
######################################################################



print(`First Written: March 31, 2023: tested for Maple 2020 `):
print(`This Version  : April  2023  `):
print():
print(`This is  Marko.txt, a Maple package`):
print(`accompanying Shalosh B. Ekhad and Doron Zeilberger's  article: `):
print(` Counting Clean Words According to The Number of Their Clean Neighbors `):
print(`IN FOND MEMORY OF MARKO PETKOVSEK (1955-2023)`):
print():
print(`The most current version is available on WWW at:`):
print(` http://sites.math.rutgers.edu/~zeilberg/tokhniot/Marko.txt .`):
print(`Please report all bugs to: DoronZeil at gmail dot com .`):

print():
print(`For general help, and a list of the MAIN functions,`):
print(` type "ezra();". For specific help type "ezra(procedure_name);" `):

print(`For a list of the supporting functions type: ezra1();`):
print(`For a list of the Story functions type: ezraSt();`):
print(`For a list of the Checking functions type: ezraC();`):
print():

with(combinat):




ezraSt:=proc()
if args=NULL then

print(`The Story procedures are: Paper1, Paper2, Paper3, Theor `):

else
ezra(args):

fi:

end:

ezraC:=proc()
if args=NULL then

print(`The Checking procedures are:  CheckWtEg,  CheckWtEs, SeqBF, SeqBFsSU, SeqBFs`):

else
ezra(args):

fi:


end:


ezra1:=proc()
if args=NULL then

print(`The SUPPORTING procedures are:  AsyRat, AvNeiNu,  CleanWords, CleanWordsSta, CleanWordsStu, ConvM, DersL, IsClean, NCN, TransTab, TransTab1, Words, WordsECg, WtG `):

else
ezra(args):

fi:


end:

ezra:=proc()
if args=NULL then

print(` Marko.txt: A Maple package for  Investigating Clean Words according to their number of clean neighbors`):

print(`The MAIN procedures are: AvNei,  Contest, StatNei, WtEg, WtEs`):
print(``):



elif nargs=1 and args[1]=AsyRat then
print(`AsyRat(R,x,a,n): inputs a rational function R of x, and a symbol for the smallest root of the denominator, outputs`):
print(`the leading asymptotics of the n-th coeff. of the MacLaurin expansion of R. Try:`):
print(`AsyRat(x/(1-x-x^2)^3,x,a,n);`):

elif nargs=1 and args[1]=AvNei then
print(`AvNei(A,M,x,MaxK,a): Asymptotic number of clean neighbors divided by (n+1), in terms of the symbol a indicating the smallest real root of the denominator, followed by the floating point. Try`):
print(`AvNei({0,1},{[1,1]},x,6,a);`):

elif nargs=1 and args[1]=AvNeiNu then
print(`AvNeiNu(A,M,MaxK,K): Appx.Asymptotic number of clean neighbors divided by (n+1), done numerically using the K-th coeff. of the Taylor expansion. Try:`):
print(`AvNeiNu({0,1},{[1,1]},6,1000);`):

elif nargs=1 and args[1]=CheckWtEg then
print(`CheckWtEg(A,M,K,K1): Checks WtEg(A,M,x,y,t,K1) for up to K terms. Try`):
print(`CheckWtEg({0,1},{[1,1]},8,4);`):


elif nargs=1 and args[1]=CheckWtEs then
print(`CheckWtEs(A,M,K,K1): Checks WtEs(A,M,y,x,K) for up to K1 terms. Try`):
print(`CheckWtEs({0,1},{[1,1]},3,16);`):


elif nargs=1 and args[1]=CleanWords then
print(`CleanWords(A,M,k): The set of words of length k in the alphabet A avoding the mistakes in M as factors. Try:`):
print(`CleanWords({0,1},{[1,1]},5);`):

elif nargs=1 and args[1]=CleanWordsSta then
print(`CleanWordsSta(A,M,k,L): The set of words of length k in the alphabet A avoding the mistakes in M as factors and starting with L. Try:`):
print(`CleanWordsSta({0,1},{[1,1]},5,[1,0]);`):

elif nargs=1 and args[1]=CleanWordsStu then
print(`CleanWordsStu(A,M,k): The set of words of length k in the alphabet A avoding the mistakes in M as factors, by really stupid brute force. Only for checking. Try:`):
print(`CleanWordsStu({0,1},{[1,1]},5);`):

elif nargs=1 and args[1]=Contest then
print(`Contest(R,k,x,y,G,MaxK,K0): Inputs pos. integers R, tries to find WtEs(A,M,y,x,MaxK) (q.v.) for ALL sets of G mistakes in the alphabet {1,2,..,R} of length k`):
print(`it also finds AvNeNu(A,M,MaxK,K0)  and ranks them accordingly. Try:`):
print(`Contest(2,2,x,y,1,6,100);`):

elif nargs=1 and args[1]=ConvM then
print(`ConvM(L): Given a list of numbers indicating [NuOfElements, SumOfRV,SecondMoment,..,K-th moment, outputs, in floating point`):
print(`[av,variance, skewness,kurtosis,...,K-th moment). Try`):
print(`ConM([5,10,100,50]);`):

elif nargs=1 and args[1]=DersL then
print(`DersL(A,M,k,L): inputs an alphabet A, a set of forbidden factors M, a prefix L and an integer k, outputs the set of`):
print(`NCN(w1,A,M)-NCN([op(2..nops(w1),w1)],A,M) over all w1 in CleanWordsSta(A,M,k,L). Try:`):
print(`DersL({0,1},{[1,1]},[1,0,1],5);`):


elif nargs=1 and args[1]=IsClean then
print(`IsClean(w,M): inputs a word w and a set of bad words M, decides whether w contains as a factor (i.e. consecutive subword) any member of M`):

elif nargs=1 and args[1]=NCN then
print(`NCN(w,A,M): Given a word w in the alphabet A and a set of mistakes M, finds the number of neighbors that are also clean. Try:`):
print(`NCN([0,1,0,1],{0,1},{[1,1]});`):


elif nargs=1 and args[1]=Paper1 then
print(`Paper1(K,x,y): The weight-enumerators of words in {0,1} avoding as a factor r repetitions of 1, according to length (given by x) and number of clean neighbors (given by y)`):
print(`and from r=2 to r=K, along with their limiting average number of clean neighbors. Try:`):
print(`Paper1(5,x,y);`):


elif nargs=1 and args[1]=Paper2 then
print(`Paper2(K,x,y): The weight-enumerators of words in {0,1} avoding as a factor r repetitions of 1, and also avoiding r repetitions of 0, according to length (given by x) and number of clean neighbors (given by y)`):
print(`and from r=2 to r=K, along with their limiting average number of clean neighbors. Try:`):
print(`Paper2(5,x,y);`):

elif nargs=1 and args[1]=Paper3 then
print(`Paper3(K,x,y): Verbose form of Contest(R,k,x,y,G,MaxK,K0) (q,v.). Try:`):
print(`Paper3(2,2,x,y,1,10,500);`):

elif nargs=1 and args[1]=SeqBF then
print(`SeqBF(A,M,x,y,K): The first K terms of the weight-enumerators of clean words. Try:`):
print(`SeqBF({0,1},{[1,1]},x,y,6);`):


elif nargs=1 and args[1]=SeqBFsSU then
print(`SeqBFsSU(A,M,y,K,SU): The first K terms of the weight-enumerators of clean words that end with SU. Try:`):
print(`SeqBFsSU({0,1},{[1,1]},y,6,[1]);`):


elif nargs=1 and args[1]=SeqBFs then
print(`SeqBFs(A,M,y,K): The first K terms of the weight-enumerators of clean words only keeping track of the number of clean neighbors. Try:`):
print(`SeqBFs({0,1},{[1,1]},y,6);`):


elif nargs=1 and args[1]=StatNei then
print(`StatNei(A,M,MaxK,K1,K2): inputs an alphabet A, a list of forbidden consecutive subwords M, a pos. integer parameter MaxK (make it big enough), a large pos. integer K1,`):
print(`and a moderate pos. integer K2 outputs the list of lists of length 10, giving for n=K1,...,n=K1+10, the (i) expectated number of clean  neighbors for a random clean word of length n/(n+1)`):
print(`(ii) the variance/(n+1) (iii) the 3rd-through the K2-th scaled moments. Try:`):
print(`StatNei({0,1},{[1,1]},4,100,4);`):


elif nargs=1 and args[1]=Theor then
print(`Theor(A,M,K,x,y): The weight-enumerators of words in the alphabet A avoding as a factor the members of M (where they all have the same length) according to length (given by x) and number of clean neighbors (given by y). Try:`):
print(`K is a guessing parameter`):
print(`Theor({0,1},{[1,0,1,0],[0,1,0,1]},6,x,y);`):

elif nargs=1 and args[1]=TransTab1 then
print(`TransTab1(A,M,k1,k2): Inputs an alphabet A, a set of mistakes M, pos. integers k1 and k2 (with k2>=k1+2), outputs`):
print(`(i) The set of clean words with of length k1, let's call it G`):
print(`(ii) If possible a table for each member w, of G, and each letter in a,`):
print(`DersL(A,M,k2,[op(w),a])[1] if  DersL(A,M,k2,[op(w),a]) is a singleton, otherwise it returns FAIL. Try:`):
print(`TransTab1({0,1},{[1,1]},2,7);`):

elif nargs=1 and args[1]=TransTab then
print(`TransTab(A,M,MaxK): Inputs an alphabet A, a set of mistakes M, pos. MaxK, tries to find the transition table for A and M. Try:`):
print(`TransTab({0,1},{[1,1]},6);`):

elif nargs=1 and args[1]=Words then
print(`Words(A,k): Inputs an alphabet A and a non-neg. integer k, outputs the set of words in A with k letters. Try:`):
print(`Words({0,1},5);`):

elif nargs=1 and args[1]=WordsEC then
print(`WordsEC(R,k): the set of words of length k in {1,..,R} but up to equivalnce under permutation. Try:`):
print(`WordsEC(3,4);`):


elif nargs=1 and args[1]=WordsECg then
print(`WordsECg(R,k): the set of words of length k in {1,..,R} but up to equivalnce under permutation. It gives one represnative. Try:`):
print(`WordsECg(3,4,1);`):

elif nargs=1 and args[1]=WtEg then
print(`WtEg(A,M,x,y,t,MaxK): Inputs an alphabet A, a set of words of the same length M, and variables x,y,t, outputs the The rational function in (x,y,t) such its MacLaurin expansion w.r.t. t, the coeff. of t^k`):
print(`is the weight enumerator of all words of length k that avoid the members of M as consecutive subwords (i.e. factors) according to the weight product of x[w[i]] times y to the power the number of clean neighbors`):
print(`(i.e. Hamming distance 1). Try:`):
print(`WtEg({0,1},{[1,1]},x,y,t,4);`):


elif nargs=1 and args[1]=WtEs then
print(`WtEs(A,M,y,x,MaxK): Inputs an alphabet A, a set of words of the same length M, and variables y,x, outputs the The rational function in (y,x) such its MacLaurin expansion w.r.t. x, the coeff. of x^k`):
print(`is the weight enumerator of all words of length k that avoid the members of M as consecutive subwords (i.e. factors) according to the weight y to the power the number of clean neighbors`):
print(`(i.e. those with Hamming distance 1). Try:`):
print(`WtEs({0,1},{[1,1]},y,x,4);`):


elif nargs=1 and args[1]=WtG then
print(`WtG(w,A,M,x,y): The general weight (keeping track of the individual letters) prod(x[w[i]]*y^NCN(w,A,M))*t^length(w). Try:`):
print(`WtG([0,1,0,1],{0,1},{[1,1]},x,y,t);`):

else
 print(`There is no such thing as`, args):

fi:


end:


#Words(A,k): Inputs an alphabet A and a non-neg. integer k, outputs the set of words in A with k letters. Try:
#Words({0,1},5);
Words:=proc(A,k) local gu,a,gu1:option remember: if k=0 then  RETURN({[]}): else gu:=Words(A,k-1): RETURN({seq(seq([a,op(gu1)],a in A),gu1 in gu)}):fi: end:

#IsClean(w,M): inputs a word w and a set of bad words M, decides whether w contains as a factor (i.e. consecutive subword) any member of M
IsClean:=proc(w,M) local i,m:
for m in M do
for i from 1 to nops(w)-nops(m)+1 do
 if [op(i..i+nops(m)-1,w)]=m then
  RETURN(false):
 fi:
od:
od:
true:
end:

#CleanWordsStu(A,M,k): The set of words of length k in the alphabet A avoding the mistakes in M as factors, by really stupid brute force. Only for checking. Try:
#CleanWordsStu({0,1},{[1,1]},5);
CleanWordsStu:=proc(A,M,k) local lu,lu1,gu:
lu:=Words(A,k):
gu:={}:
for lu1 in lu do
 if IsClean(lu1,M) then
  gu:=gu union {lu1}:
 fi:
od:
gu:
end:


#CleanWords(A,M,k): The set of words of length k in the alphabet A avoding the mistakes in M as factors. Try:
#CleanWords({0,1},{[1,1]},5);
CleanWords:=proc(A,M,k) local lu,lu1,gu,a,w1:

if k=0 then
 RETURN({[]}):
fi:
lu:=CleanWords(A,M,k-1):
gu:={}:

for lu1 in lu do
 for a in A do
   w1:=[a,op(lu1)]:
   if IsClean(w1,M) then
      gu:=gu union {w1}:
   fi:
 od:
od:

gu:

end:



#CleanWordsBF(A,M,k): The set of words of length k in the alphabet A avoding the mistakes in M as factors. Try:
#CleanWordsBF({0,1},{[1,1]},5);
CleanWordsBF:=proc(A,M,k) local lu,lu1,gu:
lu:=Words(A,k):
if k=0 then
 RETURN({[]}):
fi:
gu:={}:
for lu1 in lu do
if IsClean(lu1,M) then
 gu:=gu union {lu1}:
fi:
od:

gu:

end:



#CleanWordsSta(A,M,k,L): The set of words of length k in the alphabet A avoding the mistakes in M as factors and starting with L. Try:
#CleanWordsSta({0,1},{[1,1]},5,[1,0]);
CleanWordsSta:=proc(A,M,k,L) local lu,lu1,gu,a,w1:

if k<nops(L) then
 RETURN({}):
fi:

if not IsClean(L,M) then
  RETURN({}):
fi:

if k=nops(L) then
 RETURN({L}):
fi:

lu:=CleanWordsSta(A,M,k-1,L):
gu:={}:

for lu1 in lu do
 for a in A do
   w1:=[op(lu1),a]:
   if IsClean(w1,M) then
      gu:=gu union {w1}:
   fi:
 od:
od:

gu:



end:





#NCN(w,A,M): Given a word w in the alphabet A and a set of mistakes M, finds the number of neighbors that are also clean. Try:
#NCN([0,1,0,1],{0,1},{[1,1]});
NCN:=proc(w,A,M) local a,w1,co,i:

if not IsClean(w,M) then
 print(w, `is already bad `):
 RETURN(FAIL):
fi:

co:=0:

for i from 1 to nops(w) do
 for a in A minus {w[i]} do
   w1:=[op(1..i-1,w),a,op(i+1..nops(w),w)]:
   if IsClean(w1,M) then
      co:=co+1:
   fi:
 od:
od:
co:
end:


#WtG(w,A,M,x,y): The general weight (keeping track of the individual letters) prod(x[w[i]]*y^NCN(w,A,M))*t^length(w). Try:
#WtG([0,1,0,1],{0,1},{[1,1]},x,y,t)
WtG:=proc(w,A,M,x,y) local i:mul(x[w[i]],i=1..nops(w))*y^NCN(w,A,M):end:


#Wt(w,A,M,y): The general weight (keeping track of the individual letters) x^nops(w)*y^NCN(w,A,M))*t^length(w). Try:
#Wt([0,1,0,1],{0,1},{[1,1]},y,t)
Wt:=proc(w,A,M,y) local i:y^NCN(w,A,M):end:


#SeqBF(A,M,x,y,K): The first K terms of the weight-enumerators of clean words. Try:
#SeqBF({0,1},{[1,1]},x,y,6);
SeqBF:=proc(A,M,x,y,K) local gu,i,lu,lu1,gu1:
gu:=[]:

for i from 1 to K do
 lu:=CleanWords(A,M,i):
 gu1:=add(WtG(lu1,A,M,x,y),lu1 in lu):
 gu:=[op(gu),gu1]:
od:
gu:
end:


#SeqBFsSU(A,M,x,y,K,SU): The first K terms of the weight-enumerators of clean words that end with SU. Try:
#SeqBFsSU({0,1},{[1,1]},x,y,6,[1]);
SeqBFsSU:=proc(A,M,y,K,SU) local x,gu,i,lu,lu1,gu1,w1,a:
gu:=[]:

for i from 1 to K do
 lu:=CleanWords(A,M,i):
 gu1:=0:
  for lu1 in lu do
   w1:=[op(lu1),op(SU)]:
  if IsClean(w1,M) then
#KAKI
    gu1:=gu1+WtG(w1,A,M,x,y):
  fi:
 od:
 gu:=[op(gu),subs({seq(x[a]=1,a in A)},gu1)]:
od:
gu:
end:


#SeqBFs(A,M,y,K): The first K terms of the weight-enumerators of clean words. Try:
#SeqBF({0,1},{[1,1]},y,6);
SeqBFs:=proc(A,M,y,K) local x,a:
expand(subs({seq(x[a]=1,a in A)},SeqBF(A,M,x,y,K))):
end:


#WtEold(A,M,x,y,t): Inputs an alphabet A, a set of words of the same length M, and variables x,y,t, outputs the The rational function in (x,y,t) such its MacLaurin expansion w.r.t. t, the coeff. of t^k
#is the weight enumerator of all words of length k that avoid the members of M as consecutive subwords (i.e. factors) according to the weight product of x[w[i]] times y to the power the number of clean neighbors
#(i.e. Hamming distance 1). Try:
#WtEold({0,1},{[1,1]},x,y,t);
WtEold:=proc(A,M,x,y,t) local w,eq,var,G,k,f,eq1,a,di,w1,w2,T,var1,S,i:


if M={} then
 RETURN(FAIL):
fi:

if M<>{} and nops({seq(nops(w), w in M)})<>1 then
 RETURN(FAIL):
fi:

k:=nops(M[1]):

G:=Words(A,k) minus M:

T:=add(f[w], w in G):




for i from 0 to k-1 do
S:=Words(A,i):
T:=T+add(WtG(w,A,M,x,y), w in S)*t^i:
od:


var:={seq(f[w], w in G)}:

eq:={}:

for w in G do
eq1:=f[w]-WtG(w,A,M,x,y)*t^k:
 for a in A do
  w1:=[op(w),a]:
 w2:=[op(2..k+1,w1)]:
 if member(w2,G) then
  di:=NCN(w1,A,M)-NCN(w2,A,M):
 eq1:=eq1-x[w[1]]*t*y^di*f[w2]:
 fi:
 od:
eq:=eq union {eq1}:
od:



var1:=solve(eq,var):

if var1=NULL then
 RETURN(FAIL):
fi:

normal(subs(var1,T)):

end:


#CheckWtEg(A,M,K,K1): Checks WtEg(A,M,x,y,t,K1) for up to K terms. Try
#CheckWtEg({0,1},{[1,1]},8,4);
CheckWtEg:=proc(A,M,K,K1) local x,y,t,gu,mu,i:
gu:=WtEg(A,M,x,y,t,K):


if gu=FAIL then
 RETURN(FAIL):
fi:

gu:=taylor(gu,t=0,K1+1):
mu:=SeqBF(A,M,x,y,K1):
evalb({seq(expand(coeff(gu,t,i)-mu[i]),i=1..K1)}={0}):
end:



#CheckWtEs(A,M,K,K1): Checks WtEs(A,M,x,y,t,K) for up to K1 terms. Try
#CheckWtEs({0,1},{[1,1]},3,12);
CheckWtEs:=proc(A,M,K,K1) local x,y,t,gu,mu,i:
gu:=WtEs(A,M,y,t,K):


if gu=FAIL then
 RETURN(FAIL):
fi:

gu:=taylor(gu,t=0,K1+1):
mu:=SeqBFs(A,M,y,K1):
evalb({seq(expand(coeff(gu,t,i)-mu[i]),i=1..K1)}={0}):
end:



#DersL(A,M,k,L): inputs an alphabet A, a set of forbidden factors M, a prefix L and an integer k, outputs the set of
#NCN(w1,A,M)-NCN([op(2..nops(w1),w1)],A,M) over all w1 in CleanWordsSta(A,M,k,L). Try:
#DersL({0,1},{[1,1]},[1,0,1],5);
DersL:=proc(A,M,k,L) local gu, w:

if not IsClean(L,A,M) then
 RETURN(FAIL):
fi:

gu:=CleanWordsSta(A,M,k,L):

{seq(NCN(w,A,M)-NCN([op(2..nops(w),w)],A,M), w in gu)}:
end:

#TransTab1(A,M,k1,k2): Inputs an alphabet A, a set of mistakes M, pos. integers k1 and k2 (with k2>=k1+2), outputs
#(i) The set of clean words with of length k1, let's call it G
#(ii) If possible a table for each member w, of G, and each letter in a,
#DersL(A,M,k2,[op(w),a])[1] if  DersL(A,M,k2,[op(w),a]) is a singleton, otherwise it returns FAIL. Try:
#TransTab1({0,1},{[1,1]},2,7);
TransTab1:=proc(A,M,k1,k2) local gu,gu1,a,L,ka,T:


if not k2>=k1+2 then
 print(k2, `must be at least`, k1+2):
 RETURN(FAIL):
fi:

gu:=CleanWords(A,M,k1):

T:={}:

for gu1 in gu do
 for a in A do
 L:=[op(gu1),a]:

if IsClean(L,M) then
 ka:=DersL(A,M,k2,L):


 if nops(ka)<>1 then
   RETURN(FAIL):
 else
 T:=T union {[gu1,a,ka[1]]}:
fi:
fi:
od:
od:
  
T:

end:


#TransTab(A,M,MaxK): Inputs an alphabet A, a set of mistakes M, pos. MaxK, tries to find the transition table for A and M. Try:
#TransTab({0,1},{[1,1]},6);
TransTab:=proc(A,M,MaxK) local k,gadol,lu:
option remember:
gadol:=nops(M[1])+1:
lu:=TransTab1(A,M,1,1+gadol):
for k from 2 to MaxK while lu=FAIL do
 lu:=TransTab1(A,M,k,k+gadol):
 if TransTab1(A,M,k,k+gadol-1)<>lu then
    lu:=FAIL:
 fi:
od:
if k=MaxK+1 then
  RETURN(FAIL):
fi:
lu:
end:

#WtEg(A,M,x,y,t,MaxK): Inputs an alphabet A, a set of words of the same length M, and variables x,y,t, outputs the The rational function in (x,y,t) such its MacLaurin expansion w.r.t. t, the coeff. of t^k
#is the weight enumerator of all words of length k that avoid the members of M as consecutive subwords (i.e. factors) according to the weight product of x[w[i]] times y to the power the number of clean neighbors
#(i.e. Hamming distance 1). Try:
#WtEg({0,1},{[1,1]},x,y,t,4);
WtEg:=proc(A,M,x,y,t,MaxK) local k,lu,S,T,lu1,eq,var,f,eq1,di,a,w,TOT,S1,i,var1,gu,mu,gadol:
option remember:
lu:=TransTab(A,M,MaxK):

k:=nops(lu[1][1]):

if lu=FAIL then
 RETURN(FAIL):
fi:

S:=CleanWords(A,M,k):



for lu1 in lu do
T[lu1[1],lu1[2]]:=lu1[3]:
od:



var:={seq(f[w],w in S)}:
eq:={}:

for w in S do
 eq1:=f[w]-WtG(w,A,M,x,y)*t^nops(w):
 for a in A do
  if IsClean([op(w),a],M) then
     di:=T[w,a]:
   eq1:=eq1-t*x[w[1]]*f[[op(2..nops(w),w),a]]*y^di:
  fi:
 od:
eq:=eq union {eq1}:
od:




TOT:=0:

for i from 0 to k-1 do
 S1:=CleanWords(A,M,i):
 TOT:=TOT+add(WtG(w,A,M,x,y), w in S1)*t^i:
od:
TOT:=TOT+add(f[w], w in S):

var1:=solve(eq,var):

if var1=NULL then
 RETURN(FAIL):
fi:

TOT:=normal(subs(var1,TOT)):



mu:=SeqBF(A,M,x,y,k+4):
gu:=taylor(TOT,t=0,k+5):

if [seq(expand(coeff(gu,t,i)),i=1..k+4)]<>mu then

 RETURN(FAIL):
fi:

TOT:

end:





WtEsOld:=proc(A,M,y,t,MaxK) local gu,x,a:
gu:=WtEg(A,M,x,y,t,MaxK):
normal(subs({seq(x[a]=1,a in A)},gu)):
end:




#WtEsSlow(A,M,y,t,MaxK): Inputs an alphabet A, a set of words of the same length M, and variables x,y,t, outputs the The rational function in (x,y,t) such its MacLaurin expansion w.r.t. t, the coeff. of t^k
#is the weight enumerator of all words of length k that avoid the members of M as consecutive subwords (i.e. factors) according to the weight product of x[w[i]] times y to the power the number of clean neighbors
#(i.e. Hamming distance 1). Try:
#WtEsSlow({0,1},{[1,1]},x,y,t,4);
WtEsSlow:=proc(A,M,y,t,MaxK) local k,lu,S,T,lu1,eq,var,f,eq1,di,a,w,TOT,S1,i,var1,mu,gu,gadol:
option remember:
lu:=TransTab(A,M,1,6):
gadol:=nops(M[1])+1:


for k from 2 to MaxK while lu=FAIL do


 lu:=TransTab(A,M,k,k+gadol):


 if TransTab(A,M,k,k+gadol-1)<>lu then
    lu:=FAIL:
 fi:
od:
if k=MaxK+1 then
  RETURN(FAIL):
fi:
S:=CleanWords(A,M,k-1):


for lu1 in lu do
T[lu1[1],lu1[2]]:=lu1[3]:
od:

var:={seq(f[w],w in S)}:
eq:={}:

for w in S do
 eq1:=f[w]-Wt(w,A,M,y)*t^nops(w):
 for a in A do
  if IsClean([op(w),a],M) then
     di:=T[w,a]:
   eq1:=eq1-t*f[[op(2..nops(w),w),a]]*y^di:
  fi:
 od:
eq:=eq union {eq1}:
od:


TOT:=0:

for i from 0 to k-2 do
 S1:=CleanWords(A,M,i):
 TOT:=TOT+add(Wt(w,A,M,y), w in S1)*t^i:
od:
TOT:=TOT+add(f[w], w in S):

var1:=solve(eq,var):

if var1=NULL then
 RETURN(FAIL):
fi:

TOT:=normal(subs(var1,TOT)):



#mu:=SeqBFs(A,M,y,k+gadol+3):
mu:=SeqBFs(A,M,y,k+3):

gu:=taylor(TOT,t=0,k+4):
if [seq(expand(coeff(gu,t,i)),i=1..k+3)]<>mu then
 RETURN(FAIL):
fi:

if normal(TOT-1/(1-nops(A)*t*y))=0 then
 RETURN(FAIL):
fi:

TOT:

end:


#WtEs(A,M,y,t,MaxK): Inputs an alphabet A, a set of words of the same length M, and variables y,t, outputs the The rational function in (y,t) such its MacLaurin expansion w.r.t. t, the coeff. of t^k
#is the weight enumerator of all words of length k that avoid the members of M as consecutive subwords (i.e. factors) according to the weight  y to the power the number of clean neighbors
#(i.e. those with Hamming distance 1). Try:
#WtEs({0,1},{[1,1]},x,y,t,4);
WtEs:=proc(A,M,y,t,MaxK) local k,lu,S,T,lu1,eq,var,f,eq1,di,a,w,TOT,S1,i,var1,mu,gu,gadol:
option remember:
lu:=TransTab(A,M,MaxK):
k:=nops(lu[1][1]):
if lu=FAIL then
 RETURN(FAIL):
fi:

S:=CleanWords(A,M,k):


for lu1 in lu do
T[lu1[1],lu1[2]]:=lu1[3]:
od:

var:={seq(f[w],w in S)}:
eq:={}:

for w in S do
 eq1:=f[w]-Wt(w,A,M,y)*t^nops(w):
 for a in A do
  if IsClean([op(w),a],M) then
     di:=T[w,a]:
   eq1:=eq1-t*f[[op(2..nops(w),w),a]]*y^di:
  fi:
 od:
eq:=eq union {eq1}:
od:


TOT:=0:

for i from 0 to k-1 do
 S1:=CleanWords(A,M,i):
 TOT:=TOT+add(Wt(w,A,M,y), w in S1)*t^i:
od:
TOT:=TOT+add(f[w], w in S):

var1:=solve(eq,var):

if var1=NULL then
 RETURN(FAIL):
fi:

TOT:=normal(subs(var1,TOT)):

if normal(TOT-1/(1-nops(A)*t*y))=0 then
 RETURN(FAIL):
fi:


mu:=SeqBFs(A,M,y,k+4):
gu:=taylor(TOT,t=0,k+5):
if [seq(expand(coeff(gu,t,i)),i=1..k+4)]<>mu then
 RETURN(FAIL):
fi:

TOT:

TOT:

end:



#AsyRat(R,x,a,n): inputs a rational function R of x, and a symbol for the smallest root of the denominator, outputs
#the leading asymptotics of the n-th coeff. of the MacLaurin expansion of R. Try:
#AsyRat(x/(1-x-x^2)^3,x,a,n);
AsyRat:=proc(R,x,a,n) local R1,P,Q,k,Den:

#R1:=factor(R):
R1:=R:
P:=numer(R1):
Den:=denom(R1):




if not type(Den,`^`) then
Q:=Den:
k:=1:
else
Q:=op(1,Den):
k:=op(2,Den):
fi:

(-1)^k*subs(x=a,P)/subs(x=a,diff(Q,x))^k/a^k*expand(binomial(k+n-1,k-1))*(1/a)^n:
end:


#AvNei(A,M,x,MaxK,a): Asymptotic number of clean neighbors divided by (n+1), in terms of the symbol a indicating the smallest real root of the denominator, followed by the floating point
AvNei:=proc(A,M,x,MaxK,a) local y,f,f0,f1,lu0,lu1,lu,n,p,gu,cha,si, hal,i:

f:=WtEs(A,M,y,x,MaxK):

if f=FAIL then
 RETURN(FAIL):
fi:

f0:=subs(y=1,f):
f1:=normal(subs(y=1,diff(f,y))):
lu0:=AsyRat(f0,x,a,n):
lu1:=AsyRat(f1,x,a,n):

lu:=simplify(lu1/lu0):
lu:=normal(lu/(n+1)):

p:=denom(f0):



gu:=[solve(p,x)]:

cha:=gu[1]:
si:=evalf(abs(gu[1])):

for i from 2 to nops(gu) do
 hal:=evalf(abs(gu[i])):
 if hal<si then
   cha:=gu[i]:
   si:=hal:
 fi:
od:

[lu,p,evalf(subs(a=cha,lu))]:

end:






#ConvM(L): Given a list of numbers indicating [NuOfElements, SumOfRV,SecondMoment,..,K-th moment, outputs, in floating point
#[av,variance, skewness,kurtosis,...,K-th moment). Try
#ConV([5,10,100,50]);
ConvM:=proc(L) local K,mu,var,M,r,j:
K:=nops(L)-1:
mu:=evalf(L[2]/L[1]):
M[1]:=mu:
var:=evalf(L[3]/L[1]-mu^2):
M[2]:=var:
 
 for r from 3 to K do
   M[r]:=evalf(add((-1)^j*binomial(r,j)*(L[r-j+1]/L[1])*mu^j,j=0..r)):
   M[r]:=evalf(M[r]/var^(r/2)):
 od:

[seq(M[j],j=1..K)]:
end:

#StatNei(A,M,MaxK,K1,K2): inputs an alphabet A, a list of forbidden consecutive subwords M, a pos. integer parameter MaxK (make it big enough), a large pos. integer K1,
#and a moderate pos. integer K2 outputs the list of lists of length 10, giving for n=K1,...,n=K1+10, the (i) expectated number of clean  neighbors for a random clean word of length n/(n+1)
#(ii) the variance/(n+1) (iii) the 3rd-through the K2-th scaled moments. Try:
#StatNei({0,1},{[1,1]},6,1000,6);
StatNei:=proc(A,M,MaxK,K1,K2) local f,i,F,lu,j,x,y,Gu,gu:
f:=WtEs(A,M,y,x,MaxK):
if f=FAIL then
 print(`Make `, MaxK, `larger `):
 RETURN(FAIL):
fi:

lu:=f:
for i from 0 to K1 do
 F[i]:=normal(subs(y=1,lu)):
 lu:=normal(y*diff(lu,y)):
od:

for i from 0 to K1 do
F[i]:=taylor(F[i],x=0,K1+20):
od:


Gu:=[]:
for j from K1 to K1+10 do
 gu:=ConvM([seq(coeff(F[i],x,j),i=0..K2)]):
 gu:=[gu[1]/(j+1),gu[2]/(j+1),op(3..nops(gu),gu)]:
 Gu:=[op(Gu),gu]:
od:

Gu:
#[seq(ConvM([seq(coeff(F[i],x,j),i=0..K2)]),j=K1..K1+10)]:
#[seq(ConvM([coeff(F[1],x,j)/(j+1),coeff(F[2],x,j)/(j+1),seq(coeff(F[i],x,j),i=3..K2)]),j=K1..K1+10)]:
end:



#Paper1(K,x,y): The weight-enumerators of words in {0,1} avoding as a factor r repetitions according to length (given by x) and number of clean neighbors (given by y)
#and from r=2 to r=K, along with their limiting average number of clean neighbors. Try:
#Paper1(5,x,y);
Paper1:=proc(K) local x,y,a,r,co,A,m,n,f,lu,t0,M,C:
t0:=0:
print(``):
A:={0,1}:
print(`Counting  Words in the Alphabet {0,1}, That Avoid r consecutive Occurences of the letter 1 According to Their Number of Clean Neighbors for r from 2 to`, K):
print(``):
print(`By Shalosh B. Ekhad `):
print(``):
co:=0:

for r from 2 to K do
M:={[1$r]}:
f:=WtEs(A,M,y,x,3*K+4) :
if f<>FAIL then
co:=co+1:
print(``):
print(`--------------------------------------------------------------------------`):
print(``):
print(cat("Theorem ", co)):
print(``):
print(`Let C(m,n) be the number of words of length n in the alphabet`, A, `avoiding `, r, `consecutive occurrences of the letter 1, and having m neighbors (i.e. Hamming distance 1) that also obey this property, then`):

print(Sum(Sum(C(m,n)*x^n*y^m,n=0..infinity),m=0..infinity)=f):
print(``):
print(`and in Maple notation`):
lprint(f):
print(``):
lu:=AvNei(A,M,x,K+3,a):
print(`As the length of the word goes to infinity, the average number of clean neighbors of a random word of length n tends to n times`, lu[1]):
print(`where a is the root of the polynomial`, lu[2], `and in decimals this is`, lu[3]):
fi:
od:
print(``):
print(`----------------------------------`):
print(``):
print(`This ends this paper that took`, time()-t0, ` seconds to produce `):


end:



#Paper2(K,x,y): The weight-enumerators of words in {0,1} avoding as a factor r repetitions of the same letter, according to length (given by x) and number of clean neighbors (given by y)
#and from r=2 to r=K, along with their limiting average number of clean neighbors. Try:
#Paper2(5,x,y);
Paper2:=proc(K) local x,y,a,r,co,A,m,n,f,lu,t0,M,C:
t0:=0:
print(``):
A:={0,1}:
print(`Counting  Words in the Alphabet {0,1}, That Avoid r consecutive Occurences of the same letter, According to Their Number of Clean Neighbors for r from 3 to`, K):
print(``):
print(`By Shalosh B. Ekhad `):
print(``):
co:=0:

for r from 3 to K do
M:={[0$r],[1$r]}:
f:=WtEs(A,M,y,x,K+3) :
if f<>FAIL then
co:=co+1:
print(``):
print(cat("Theorem ", co)):
print(``):
print(`Let C(m,n) be the number of words of length n in the alphabet`, A, `avoiding strings of`, r, `consecutive occurrences of the same letter (i.e. you are not allowed`):
print(`to have consecutive substrings in the set`, M, `and having m neighbors (i.e. Hamming distance 1) that also obey this property, then`):

print(Sum(Sum(C(m,n)*x^n*y^m,n=0..infinity),m=0..infinity)=f):
print(``):
print(`and in Maple notation`):
lprint(f):
print(``):
lu:=AvNei(A,M,x,K+3,a):
print(`As the length of the word goes to infinity, the average number of clean neighbors of a random word of length n tends to n times`, lu[1]):
print(`where a is the root of the polynomial`, lu[2], `and in decimals this is`, lu[3]):
fi:
od:
print(``):
print(`----------------------------------`):
print(``):
print(`This ends this paper that took`, time()-t0, `seconds to produce `):


end:




#Theor(A,M,K,x,y): The weight-enumerators of words in the alphabet A avoding as a factor the members of M (where they all have the same length) according to length (given by x) and number of clean neighbors (given by y). Try:
#K is a guessing parameter
#Theor({0,1},{[1,0,1,0],[0,1,0,1]},6,x,y);
Theor:=proc(A,M,K,x,y) local C,m,n,a,f,lu,t0:

t0:=0:

f:=WtEs(A,M,y,x,K) :
if f=FAIL then
 print(`Try to make`, K, `bigger `):
 RETURN(FAIL):
fi:

print(``):
print(`Counting  Words in the Alphabet`, A, `That Avoid as consecutive subwords members of the set of words`, M):
print(``):
print(`By Shalosh B. Ekhad `):
print(``):


print(``):
print(`Theorem: Let C(m,n) be the number of words of length n in the alphabet`, A, `avoiding consecutive subsords that belong to`,M):
print(`and having m neighbors (i.e. Hamming distance 1) that also obey this condition, then`):

print(Sum(Sum(C(m,n)*x^n*y^m,n=0..infinity),m=0..infinity)=f):
print(``):
print(`and in Maple notation`):
lprint(f):
print(``):
lu:=AvNei(A,M,x,K,a):
print(`As the length of the word goes to infinity, the average number of clean neighbors of a random word of length n tends to n times`, lu[1]):
print(`where a is the smallest root (in absolute value) of the polynomial`, lu[2], `and in decimals this is`, coeff(lu[3],I,0)):


print(``):
print(`----------------------------------`):
print(``):
print(`This ends this paper that took`, time()-t0, `seconds to produce `):


end:


#WordsEC(R,k): the set of words of length k in {1,..,R} but up to equivalence under permutation. Try:
#WordsEC(3,4);
WordsEC:=proc(R,k) local i,lu,lu1,gu,P,temu,pi:
lu:=Words({seq(i,i=1..R)},k) :
P:=permute(R):
gu:={}:

while lu<>{} do
lu1:=lu[1]:
temu:={seq(subs({seq(i=pi[i],i=1..R)},lu1),pi in P)}:
gu:=gu union {temu}:
lu:=lu minus temu:
od:
gu:

end:


#WordsECg(R,k,G): the set of set of words of cardinality G of {1,..,R} but up to equivalence under permutation. It gives one representative. Try:
#WordsECg(3,4,2);
WordsECg:=proc(R,k,G) local i,lu,lu1,gu,P,temu,pi:
lu:=choose(Words({seq(i,i=1..R)},k),G) :

P:=permute(R):
gu:={}:

while lu<>{} do
lu1:=lu[1]:
temu:={seq(subs({seq(i=pi[i],i=1..R)},lu1),pi in P)}:
gu:=gu union {temu}:
lu:=lu minus temu:
od:

{seq(gu[i][1],i=1..nops(gu))}:

end:





#AvNeiNu(A,M,MaxK,K): Appx.Asymptotic number of clean neighbors divided by (n+1), done numerically using the K-th coeff. of the Taylor expansion. Try:
#AvNeiNu({0,1},{[1,1]},6,1000);
AvNeiNu:=proc(A,M,MaxK,K) local f,x,y,f0,f1:

f:=WtEs(A,M,y,x,MaxK) :
if f=FAIL then
 RETURN(FAIL):
fi:

f0:=normal(subs(y=1,f)):

f1:=normal(subs(y=1,diff(f,y))):
evalf(coeff(taylor(f1,x=0,K+1),x,K)/coeff(taylor(f0,x=0,K+1),x,K)/(K+1)):

end:

#Contest(R,k,x,y,G,MaxK,K0): Inputs pos. integers R, tries to find WtEs(A,M,y,x,MaxK) (q.v.) for ALL sets of G mistakes in the alphabet {1,2,..,R} of length k
#it also finds AvNeNu(A,M,MaxK,K0) and ranks them accordingly. Try:
#Contest(2,2,x,y,,1,10,500);
Contest:=proc(R,k,x,y,G,MaxK,K0) local A,gu,gu1,T1,T2,f,NG,i,f0,f1,T2inv,GOO,g,GOOV,x1,a1,i1:
gu:=WordsECg(R,k,G):
A:={seq(i,i=1..R)}:

GOO:={}:
NG:={}:
for gu1 in gu do
 f:=WtEs(A,gu1,y,x,MaxK):
 if f=FAIL then
  NG:=NG union {gu1}:
else
 GOO:=GOO union {gu1}:
 T1[gu1]:=f:

f0:=normal(subs(y=1,f)):

f1:=normal(subs(y=1,diff(f,y))):
T2[gu1]:=evalf(coeff(taylor(f1,x=0,K0+1),x,K0)/coeff(taylor(f0,x=0,K0+1),x,K0)/(K0+1)):

 fi:
od:


GOOV:={seq(T2[g],g in GOO)}:

for x1 in GOOV do
T2inv[x1]:={}:
od:

for gu1 in GOO do
 T2inv[T2[gu1]]:=T2inv[T2[gu1]] union {gu1}:
od:


op(T1),op(T2),NG,GOOV,op(T2inv):

GOOV:=sort(convert(GOOV,list),`>`):

GOO:=[seq(op(T2inv[GOOV[a1]]),a1=1..nops(GOOV))]:

[seq([GOO[i1],T1[GOO[i1]],T2[GOO[i1]]],i1=1..nops(GOO))],NG:
end:






#Paper3(K,x,y): Verbose form of Contest(R,k,x,y,G,MaxK,K0) (q,v.). Try:
#Paper3(2,2,x,y,,1,10,500);
Paper3:=proc(R,k,x,y,G,MaxK,K0) local gu,i,A,a,t0,lu,m,n:
gu:=Contest(R,k,x,y,G,MaxK,K0)[1]:

t0:=0:
print(``):
A:={seq(i,i=1..R)}:

print(`Counting  Words in the Alphabet`,A,  `That Avoid A Certain set of `, G, `consecutive subwords of length`, k, `For all the possibilities, According to their Number of Clean Neighbors`):
print(``):
print(`By Shalosh B. Ekhad `):
print(``):

for i from 1 to nops(gu) do
print(``):
print(`------------------------------------------------`):
print(``):
print(cat("Theorem Number ", i)):
print(``):

print(`Let C(m,n) be the number of words of length n in the alphabet`, A, `avoiding consecutive substrings in the set`, gu[i][1], `and having m neighbors (i.e. Hamming distance 1) that also obey this property, then`):
print(Sum(Sum(C(m,n)*x^n*y^m,n=0..infinity),m=0..infinity)=gu[i][2]):
print(``):
print(`and in Maple notation`):
lprint(gu[i][2]):
if gu[i][3]>0.1 then
lu:=AvNei(A,gu[i][1],x,MaxK,a):
if type(lu[3],float) then
 print(`As the length of the word goes to infinity, the average number of clean neighbors of a random word of length n tends to n times`, lu[1]):
 print(`where a is the root of the polynomial`, lu[2], `and in decimals this is`, lu[3]):
print(``):
 print(`BTW the ratio for words with`, K0, `letters is`, gu[i][3]):
else
print(``):
 print(`As the length of the word goes to infinity, the average number of clean neighbors of a random word of length n tends to n times, approximaley to`, gu[i][3]):
print(`[This estimate was obtained by looking at words of length`, K0, ` ] `):
fi:
fi:
od:
print(``):
print(`----------------------------------`):
print(``):
print(`This ends this paper that took`, time()-t0, `seconds to produce `):


end: