######################################################################
##Compositions.txt: Save this file as  Compositions.txt              #
## To use it, stay in the                                            #
##same directory, get into Maple (by typing: maple <Enter> )         #
##and then type:  read `Compositions.txt`<Enter>                     #
##Then follow the instructions given there                           #
##                                                                   #
##Written by Doron Zeilberger, Rutgers University ,                  #
#DoronZeil at gmail dot com                                          #
######################################################################
 
#Created: Jan. 2, 2019
 
print(`Created:  Jan. 2, 2019`):
print(` This is Compositions.txt `):
print(`It is one of the Maple packages that accompany the article `):
print(`  The "Monkey Typing Shakespeare" Problem for Compositions`):
print(`by Shalosh B. Ekhad and Doron Zeilberger`):
print(` available from Zeilberger's website`):
print(``):
print(`Please report bugs to DoronZeil at gmail dot com `):
print(``):
 print(`The most current version of this  package and paper`):
 print(` are  available from`):
 print(`http://sites.math.rutgers.edu/~zeilberg/mamarim/mamarimhtml/comp.html  .`):


print(`---------------------------------------`):
 print(`For a list of the Supporting procedures type ezra1();, for help with`):
 print(`a specific procedure, type ezra(procedure_name);   .`):
 print(``):
print(`---------------------------------------`):

print(`---------------------------------------`):
 print(`For a list of the Multi-Composition procedures type ezraM();, for help with`):
 print(`a specific procedure, type ezra(procedure_name);   .`):
 print(``):
print(`---------------------------------------`):


print(`---------------------------------------`):
 print(`For a list of the Cluster procedures type ezraC();, for help with`):
 print(`a specific procedure, type ezra(procedure_name);   .`):
 print(``):
print(`---------------------------------------`):


print(`---------------------------------------`):
 print(`For a list of the MAIN procedures type ezra();, for help with`):
 print(`a specific procedure, type ezra(procedure_name);   .`):
 print(``):
print(`---------------------------------------`):

with(combinat):

ezraM:=proc()

if args=NULL then
 print(` The multi-compositions procedures are:  AbaSet, AvotSet, GFset, GoodCseqD, InitialStatesS, InitialWt, InitialWtSt, Reduce1, Reduce11`):
 print(`StatesCset  `):

else
ezra(args):
fi:

end:

ezraC:=proc()

if args=NULL then
 print(` The the cluster procedures are: Aba, Avot, Clus, Clus1, MaxC, StatesC `):
 print(``):

else
ezra(args):
fi:

end:

ezra1:=proc()

if args=NULL then
 print(` The supporting procedures are: Av,  CheckGFone, CheckGFset, CheckGFset1, CheckGFset11, Comps, Comps1, IsUni, `):
 print(` GFone, GFoneU, GoodC, GoodC1, GoodCnaive, GoodC1naive, GoodCseqD,  IsCon, IsLess1, `):
 print(` Ove, Ove1, Skyline, Taur, Wt `):
 print(``):

else
ezra(args):
fi:

end:

ezra:=proc()

if args=NULL then
 print(`The main procedures are:   Info,  InfoV, GFset, SuperInfo1, SuperInfo1t, SuperInfo1tV, SuperInfo1V,  SuperInfo2kr, SuperInfo2krV   `):
 print(` `):


elif nops([args])=1 and op(1,[args])=Aba then
print(`Aba(C,S,r,x,t): inputs a composition C, a state S, and r between 2 and nops(C), outputs the new state and the`):
print(`change in weight in terms of x and t. Try: `):
print(`Aba([2,1,2,1,2],[2,1,2,2],4,x,t);`):


elif nops([args])=1 and op(1,[args])=AbaSet then
print(`AbaSet(SetC1,S,r,x,t): inputs a set of compositions SetC1, a state S, and r between 2 and nops(Skyline(SetC)), outputs the new state and the`):
print(`change in weight in terms of x and t. Try: `):
print(`AbaSet({[2,1,2,1,2]},[2,1,2,2],4,x,t);`):

elif nops([args])=1 and op(1,[args])=Av then
print(`Av(C,SetC): does the composition C avoid the set of compositions SetC. Try:`):
print(`Av([3,4,2],{[2,5],[4,1]});`):

elif nops([args])=1 and op(1,[args])=Avot then
print(`Avot(C,S): all the fathers of the state S with composition C. Try:`):
print(`Avot([2,1,2,1,2],[2,1,2,2]);`):

elif nops([args])=1 and op(1,[args])=AvotSet then
print(`AvotSet(SetC,S): all the fathers of the state S with a set of compositions SetC. Try:`):
print(`AvotSet({[2,1,2,1,2]},[2,1,2,2]);`):

elif nops([args])=1 and op(1,[args])=CheckGFone then
print(`CheckGFone(C,N): checks GFone(C,x) up to N terms. Try:`):
print(`CheckGFone([2,1,2],12);`):

elif nops([args])=1 and op(1,[args])=CheckGFset then
print(`CheckGFset(R,N): checks GFset(SetC,x) up to N terms for all subsets SetC of the set of compositions of <=r with the same number of  parts. Try:`):
print(`CheckGFset(3,10);`):

elif nops([args])=1 and op(1,[args])=CheckGFset1 then
print(`CheckGFset1(r,k,N): checks GFset(SetC,x) up to N terms for all subsets SetC of the set of compositions of r with k parts. Try:`):
print(`CheckGFset1(3,2,10);`):

elif nops([args])=1 and op(1,[args])=CheckGFset11 then
print(`CheckGFset11(SetC,N): checks GFset(SetC,x) up to N terms. Try:`):
print(`CheckGFset11({[2,1,2]},12);`):

elif nops([args])=1 and op(1,[args])=Clus then
print(`Clus(C,k): The set of clusters with one composition C of height k. Try: `):
print(`Clus([2,1,2],3);`):

elif nops([args])=1 and op(1,[args])=Clus1 then
print(`Clus1(C,k,m): The set of k by m clusters with one composition C. Try: `):
print(`Clus1([2,1,2],3,6);`):



elif nops([args])=1 and op(1,[args])=Comps then
print(`Comps(n,k): the set of compositions of n Try:`):
print(`Comps(5);`):

elif nops([args])=1 and op(1,[args])=Comps1 then
print(`Comps1(n,k): the set of compositions of n into exactly k parts. Try:`):
print(`Comps1(5,2);`):

elif nops([args])=1 and op(1,[args])=GFset then
print(`GFset(S,x): inputs a set of compositions S and outputs the generating function enumerating compositions that avoid`):
print(`the members of S as sub-compositions. Try:`):
print(`GFset({[2,3,2],[4,1,4]},x);`):

elif nops([args])=1 and op(1,[args])=GoodC1 then
print(`GoodC1(n,k,SetC): The set of compositions of n into k partsavooding the set of compositions SetC as subcompositions. Try:`):
print(`GoodC1(6,3,{[1,2,1]});`):

elif nops([args])=1 and op(1,[args])=GoodC then
print(`GoodC(n,SetC): The set of compositions of n avooding the set of compositions SetC as subcompositions.Try:`):
print(`GoodC(6,{[1,2,1]});`):

elif nops([args])=1 and op(1,[args])=GoodCnaive then
print(`GoodCnaive(n,SetC): The set of compositions of n avooding the set of compositions SetC as subcompositions. Naive approach. Try:`):
print(`GoodCnaive(6,{[1,2,1]});`):

elif nops([args])=1 and op(1,[args])=GoodC1naive then
print(`GoodC1naive(n,k,SetC): The set of compositions of n into k partsavooding the set of compositions SetC as subcompositions. Naive approach. Try:`):
print(`GoodC1naive(6,3,{[1,2,1]});`):

elif nops([args])=1 and op(1,[args])=GFone then
print(`GFone(C,x): the generating function of the sequence enumerating compositions of n avoiding the sub-composition C. Try`):
print(`GFone([2,1,2],x);`):

elif nops([args])=1 and op(1,[args])=GFoneU then
print(`GFoneU(C,x): inputs a UNIMODAL composition  C, and a variable name, x, and outputs the generating function`):
print(`of the enumerating sequence of compositions of n avooding C as a sub-composition. Try:`):
print(`GFoneU([2,3,2],x);`):

elif nops([args])=1 and op(1,[args])=GoodCseqD then
print(`GoodCseqD(N,SetC): the first N+1 terms, starting at n=0 and ending with n=N,`):
print(` of the enumerating sequence for compositions of n avoiding as sub-compositions the`):
print(`members of the set SetC, done DIRECTLY (by counting). For checking only! Try:`):
print(`GoodCseqD(10,{[2,3,2]});`):


elif nops([args])=1 and op(1,[args])=Info then
print(`Info(SetC,x,N): inputs a set of compositions of the same length, SetC, and a variable x, and a positive integer N, returns`):
print(`(i) The input set , SetC`):
print(`(ii) The generating function of the seqeunce enumerating compositions that avoid the compositions in SetC. let's call it a(n)`):
print(`(iii): The first N terms of the enumerating sequence.`):
print(`(iv) the limit of a(n+1)/a(n), as n goes to infinity, let's call it alpha`):
print(`(v) The limit of a(n)/alpha^n as n goes to infinity. Try:`):
print(`Info({[2,3,2]},x,30);`):

elif nops([args])=1 and op(1,[args])=InfoV then
print(`InfoV(SetC,x,N): verbose form of Info(SetC,x,N) (q.v.). Try: `):
print(`InfoV({[2,3,2]},x,30);`):

elif nops([args])=1 and op(1,[args])=InitialStatesS then
print(`InitialStatesS(SetC): the set of initial states for a set of compositions SetC. Try:`):
print(`InitialStatesS({[1,2,1],[4,1]});`):

elif nops([args])=1 and op(1,[args])=InitialWt then
print(`InitialWt(SetC1,x,t): The  initial weight of the set SetC1. Try:`):
print(`InitialWt({[1,3,1],[3,1,3]},x,t);`):

elif nops([args])=1 and op(1,[args])=InitialWtSt then
print(`InitialWtSt(SetC,S,x,t): The  initial weight of the state S with set of compositions SetC. Try:`):
print(`InitialWtSt({[2,1,2]},[2,1,2],x,t);`):

elif nops([args])=1 and op(1,[args])=IsCon then
print(`IsCon(C1,C2): is composition C1 contained in composition C2? Try: `):
print(` IsCon([1,2,1],[3,4,3,2]); `):

elif nops([args])=1 and op(1,[args])=IsLess1 then
print(`IsLess1(L1,L2) is the list L1 pointwise less than L2. try:`):
print(`IsLess1([3,4],[4,5]);`):

elif nops([args])=1 and op(1,[args])=IsUni then
print(`IsUni(C): Is the composition C unimodal? Try:`):
print(`IsUni([1,2,1]); IsUni([2,1,2]);`):


elif nops([args])=1 and op(1,[args])=MaxC then
print(`MaxC(M): given a rectangular matrix, outputs the row vector that is the max of each column. Try:`):
print(`MaxC([[2,1,2],[1,2,1]]);`):

elif nops([args])=1 and op(1,[args])=Ove then
print(`Ove(S1,S2,x,t): the overlaps between sets of compositions S1 and S2 try:`):
print(`Ove({[2,3,2]},{[4,3,2]},x,t);`):

elif nops([args])=1 and op(1,[args])=Ove1 then
print(`Ove1(C,x,t): the self overlap of a composition C . Try:`):
print(`Ove1([2,3,4,3,2],x,t);`):


elif nops([args])=1 and op(1,[args])=Reduce1 then
print(`Reduce1(Cset): given a set of compositions finds a minimal subset. Try:`):
print(`Reduce1({[1,2,1],[3,3,3],[6,4,3]});`):

elif nops([args])=1 and op(1,[args])=Reduce11 then
print(`Reduce11(Cset): given a set of compositions performs one step in kicking out superflous compositions. Try:`):
print(`Reduce11({[1,2,1],[3,3,3]});`):

elif nops([args])=1 and op(1,[args])=Skyline then
print(`Skyline(S): inputs a set of compositions, outputs the smallest composition that  includes them all. Try: `):
print(` Skyline({[1,2,1],[2,1,1],[4]}); `):

elif nops([args])=1 and op(1,[args])=StatesC then
print(`StatesC(C): The set of states for the clusters of the composition C. Try:`):
print(`StatesC([2,1,2]);`):

elif nops([args])=1 and op(1,[args])=StatesCset then
print(`StatesCset(SetC): The set of states for the clusters of the set of compositions SetC. Try:`):
print(`StatesCset({[2,1,2]});`):

elif nops([args])=1 and op(1,[args])=SuperInfo1 then
print(`SuperInfo1(k,x,N): Inputs a positive integer k, variable x, and positive integer N, gathers Info({C},x,N) for all`):
print(`compositions C of k. Try:`):
print(`SuperInfo1(5,x,30);`):

elif nops([args])=1 and op(1,[args])=SuperInfo1t then
print(`SuperInfo1t(k): Terse version of SuperInfo(k,x,N). Only gives the ascending list [comps, growh_constant]. Try: `):
print(`SuperInfo1t(5);`):

elif nops([args])=1 and op(1,[args])=SuperInfo1tV then
print(`SuperInfo1tV(K): Verbose form  of SuperInfo1t(k) (q.v.)  for k from 2 to K: Try:`):
print(`SuperInfo1tV(5);`):

elif nops([args])=1 and op(1,[args])=SuperInfo1V then
print(`SuperInfo1V(k,x,N): Verbose form of SuperInfo1(k,x,N) (q.v.). Try: `):
print(`SuperInfo1V(5,x,30);`):

elif nops([args])=1 and op(1,[args])=SuperInfo2kr then
print(`SuperInfo2kr(k,r,x,N): Inputs  positive integers k and r, variable x, and positive integer N, gathers Info(S,x,N) `):
print(`for all subsets of the set of compositions of k with r parts. Try:`):
print(`SuperInfo2kr(5,2,x,30);`):

elif nops([args])=1 and op(1,[args])=SuperInfo2krV then
print(`SuperInfo2krV(k,r,x,N): verbose form of SuperInfo2kr(r,k,x,N) (q.v.). Try: `):
print(`SuperInfo2krV(5,2,x,30);`):


elif nops([args])=1 and op(1,[args])=Taur then
print(`Taur(r): the composition r(r+1)...(2*r-2)(2*r-1)...(r+1)r. For example, try:`):
print(`Taur(4);`):

elif nops([args])=1 and op(1,[args])=Wt then
print(`Wt(S,x,t): the weight of a set of compositions S. Try:`):
print(`Wt({[1,2,1],[2,1,1],[4]},x,t);`):


else
print(`There is no ezra for`,args):
fi:
 
end:


#Comps1(n,k): the set of compositions of n into exactly k parts. Try:
#Comps1(5,2);

Comps1:=proc(n,k) local gu,mu1,mu,i:
option remember:
if k=0 then
 if n=0 then
  RETURN({[]}):
 else
  RETURN({}):
 fi:
fi:

if  n<k then
 RETURN({}):
fi:

if n=k then
 RETURN({[1$k]}):
fi:

if k=1 then
 RETURN({[n]}):
fi:

gu:={}:

for i from 1 to n do
 mu:=Comps1(n-i,k-1):
 gu:=gu union {seq([i,op(mu1)],mu1 in mu)}:
od:

gu:

end:


#Comps(n): the set of compositions of n. Try:
#Comps(5);
Comps:=proc(n) local k:
option remember:
if n=0 then
 RETURN({[]}):
fi:

{seq(op(Comps1(n,k)),k=1..n)}:
end:


#IsLess1(L1,L2) is the list L1 pointwise less than L2. try:
#IsLess1([3,4],[4,5]);
IsLess1:=proc(L1,L2) local i:
if nops(L1)<>nops(L2) then
 RETURN(FAIL):
fi:

for i from 1 to nops(L1) do
 if L1[i]>L2[i] then
  RETURN(false):
 fi:
od:
true:
end:


#IsCon(C1,C2): is composition C1 contained in composition C2? Try
#IsCon([1,2,1],[3,4,3,2]);
IsCon:=proc(C1,C2) local i:

if nops(C1)>nops(C2) then
 RETURN(false):
fi:

for i from 1 to nops(C2)-nops(C1)+1 do
if IsLess1(C1,[op(i..i+nops(C1)-1,C2)]) then
 RETURN(true):
fi:
od:

false:

end:



#Av(C,SetC): does the composition C avoid the set of compositions SetC. Try
#Av([3,4,2],{[2,5],[4,1]});
Av:=proc(C,SetC) local C1:

for C1 in SetC do
 if IsCon(C1,C) then
   RETURN(false):
 fi:
od:

true:

end:

#GoodCnaive(n,SetC): The set of compositions of n avooding the set of compositions SetC as subcompositions. Naive approach. Try:
#GoodCnaive(6,{[1,2,1]});

GoodCnaive:=proc(n,SetC) local mu,mu1,gu:
option remember:
mu:=Comps(n):
gu:={}:

 for mu1 in mu do
 if Av(mu1,SetC) then
   gu:=gu union {mu1}:
 fi:
od:
gu:
end:


#GoodC1naive(n,k,SetC): The set of compositions of n into k partsavooding the set of compositions SetC as subcompositions. Naive approach. Try:
#GoodC1naive(6,3,{[1,2,1]});

GoodC1naive:=proc(n,k,SetC) local mu,mu1,gu:
option remember:
mu:=Comps1(n,k):
gu:={}:

 for mu1 in mu do
 if Av(mu1,SetC) then
   gu:=gu union {mu1}:
 fi:
od:
gu:
end:


#GoodC1(n,k,SetC): the set of compositions of n into exactly k parts avoiding the set of compositions SetC . Try:
#GoodC1(5,2,{[1,2,1]});

GoodC1:=proc(n,k,SetC) local gu,mu1,mu,i,cand:
option remember:
if k=0 then
 if n=0 then
  RETURN({[]}):
 else
  RETURN({}):
 fi:
fi:

if  (n<k or n=k or k=1) then
 RETURN(GoodC1naive(n,k,SetC)):
fi:

gu:={}:

for i from 1 to n do
 mu:=GoodC1(n-i,k-1,SetC):

  for mu1 in mu do
   cand:=[i,op(mu1)]:
    if Av(cand,SetC)  then
     gu:=gu union {cand}:
    fi:
  od:
 
od:

gu:

end:

#GoodC(n,SetC): The set of compositions of n avooding the set of compositions SetC as subcompositions. Try:
#GoodC(6,{[1,2,1]});
GoodC:=proc(n,SetC) local k:
option remember:
{seq(op(GoodC1(n,k,SetC)),k=1..n)}:
end:

#GoodCseqD(N,SetC): the first N+1 terms, starting at n=0 and ending with n=N,
# of the enumerating sequence for compositions of n avoiding as sub-compositions the
#members of the set SetC, done DIRECTLY (by counting). For checking only! Try:
#GoodCseqD(10,{[2,3,2]});
GoodCseqD:=proc(N,SetC) local n1:
[1,seq(nops(GoodC(n1,SetC)),n1=1..N)]:
end:

#GFoneU(C,x inputs a variable name, t, and a composition C, and outputs the generating function
#of the enumerating sequence of compositions of n avooding C as a sub-composition. Try:
#GFoneU([2,3,2],x
GFoneU:=proc(C,x) local L,t:
option remember:

if C=[] then
 RETURN(0):
fi:

if not IsUni(C) then
 print(C, `is not unimodal `):
 RETURN(FAIL):
fi:

L:=normal(Wt({C},x,t) /(1+Ove1(C,x,t))):

L:=subs(t=1/(1-x),L):
normal(1/(1-x/(1-x)-L)):

end:


#Ove1(C,x,t): the self overlap of a composition C . Try
#Ove1([2,3,4,3,2],x,t);
Ove1:=proc(C,x,t) local k,j,i:
k:=nops(C):

if not (type(C,list) and type(x,symbol) and type(t,symbol) ) then
 print(`Bad input`):
 RETURN(FAIL):
fi:

add( t^(i-1)*x^(add(C[j],j=1..i-1)+add( max(C[j]-C[j-i+1],0)  , j=i..k)),i=2..k):
end:



#Taur(r): the composition r(r+1)...(2*r-2)(2*r-1)...(r+1)r. For example, try:
#Taur(4);
Taur:=proc(r) local i,gu:
gu:=[seq(i,i=r..2*r-2)]:
[op(gu),2*r-1,seq(gu[nops(gu)+1-i],i=1..nops(gu))]:
end:



#Skyline(S): inputs a set of compositions, outputs the smallest composition that  includes them all. Try
#Skyline({[1,2,1],[2,1,1],[4]});
Skyline:=proc(S) local s,k,S1,i:

k:=max(seq(nops(s), s in S)):

S1:={seq([op(s),0$(k-nops(s))],s in S)}:
[seq(max( seq(s[i],s in S1)),i=1..k)]:
end:


#Wt(S,x,t): the weight of a set of compositions S. Try:
#Wt({[1,2,1],[2,1,1],[4]},x,t);
Wt:=proc(S,x,t) local gu:
gu:=Skyline(S):

(-1)^nops(S)*t^nops(gu)*x^convert(gu,`+`):
end:





#Ove(S1,S2,x,t): the overlaps between sets of compositions S1 and S2 try:
#Ove({[2,3,2]},{[4,3,2]},x,t);
Ove:=proc(S1,S2,x,t) local C1,C2,k,i,j,lu:
C1:=Skyline(S1):
C2:=Skyline(S2):
k:=nops(C1):

lu:=0:

for i from 2 to k do

 if k-i+1<=nops(C2) then

  lu:=lu+t^(i-1)*x^(add(C1[j],j=1..i-1)+add( max(C1[j]-C2[j-i+1],0)  , j=i..k) ):
 else
   lu:=lu+t^(i-1)*x^(add(C1[j],j=1..i-1)+add( max(C1[j]-C2[j-i+1],0)  , j=i..k-i+1)+add( C1[j]  , j=k-i+2..k))
 fi:

od:
end:





#IsUni(C): Is the composition C unimodal? Try:
#IsUni([1,2,1]); IsUni([2,1,2]);
IsUni:=proc(C) local i:
if nops(C)<=1 then
 RETURN(true):
fi:
if C[2]<C[1] then
   RETURN(false):
fi:

for i from 2 to nops(C)-1 do
 if C[i-1]>C[i] and C[i]<C[i+1] then
   RETURN(false):
 fi:
od:
true:
end:

#Clus1(C,k,m): The set of k by m clusters with one composition. Try
#Clus1([2,1,2],3,6);
Clus1:=proc(C,k,m) local gu,r,gu1,gu11,i1:
option remember:

if not (type(C,list) and type(k, integer) and type(m,integer) and k>0 and m>0) then
 RETURN({}):
fi:

if k=1 then
 if m=nops(C) then
  RETURN({[C]}):
  else
  RETURN({}):
 fi:
fi:

gu:={}:

for r from 1 to nops(C)-1 do
 gu1:=Clus1(C,k-1,m-r):
 gu1:={seq([[op(C),0$(m-nops(C))], seq([0$(r),op(gu11[i1])],i1=1..k-1)],gu11 in gu1)}:
 gu:=gu union gu1:
od:
gu:


end:

#Clus(C,k): The set of k by m clusters with one composition with ANY m. Try
#Clus([2,1,2],3);
Clus:=proc(C,k) local m,m1,gu:
option remember:

for m from 1 while Clus1(C,k,m)={} do od:
m1:=m:

gu:={}:

for m from m1 while Clus1(C,k,m)<>{} do
 gu:=gu union Clus1(C,k,m):
od:
gu:

end:


#MaxC(M): given a rectangular matrix, outputs the row vector that is the max of each column. Try
#MaxC([[2,1,2],[1,2,1]]);
MaxC:=proc(M) local i,i1:

[seq(max(seq(M[i][i1],i=1..nops(M))),i1=1..nops(M[1]))]:

end:


#Aba(C,S,r,x,t): inputs a composition C, a state S, and r between 2 and nops(C), outputs the new state and the
#change in weight in terms of x and t. Try
#Aba([2,1,2,1,2],[2,1,2,2],4,x,t);
Aba:=proc(C,S,r,x,t) local gu,mu,i:
option remember:
if (nops(C)<>nops(S)  or r<=1 or r>nops(C)) then
 print(`Bad input`):
 RETURN(FAIL):
fi:

gu:=[[op(C),0$(r-1)],[0$(r-1),op(S)]]:
mu:=[seq(max(gu[1][i],gu[2][i]),i=1..nops(C)+r-1)]:

[[op(1..nops(S),mu)],-t^(r-1)*x^(convert(mu,`+`)-convert(S,`+`))]:
end:

#Avot(C,S): all the fathers of the state S with composition C. Try:
#Avot([2,1,2,1,2],[2,1,2,2]);
Avot:=proc(C,S) local r,x,t:
{seq(Aba(C,S,r,x,t)[1],r=2..nops(C))}:
end:



#StatesC(C): The set of states for the clusters of the composition C. Try:
#StatesC([2,1,2]);
StatesC:=proc(C) local S,gu,gu1,gu2, gua:

S:=C:

gu:={S}:

gu1:=Avot(C,S) union gu:

while gu<>gu1 do
 gu2:=gu1 union {seq(op(Avot(C,gua)),gua in gu1)}:
 gu:=gu1:
 gu1:=gu2:
od:
gu1:
end:




#GFone(C,x): the generating function of the sequence enumerating compositions of n avoiding the sub-composition C. Try
#GFone([2,1,2],x);
GFone:=proc(C,x) local gu,eq,var, T,S,gu1,gu2,r,lu,t,X,eq1,L:
option remember:
gu:=StatesC(C):

S:=C:

for  gu1 in gu do
 for gu2 in gu do
  T[gu1,gu2]:=0:
 od:
od:

for gu1 in gu do
 for r from 2 to nops(C) do
  lu:=Aba(C,gu1,r,x,t):
   T[lu[1],gu1]:=T[lu[1],gu1]+lu[2]:
 od:
od:




var:={seq(X[gu1],gu1 in gu)}:

eq1:=X[S]+x^convert(C,`+`)*t^nops(C)-add(T[S,gu1]*X[gu1],gu1 in gu):
eq:={eq1}:


for gu1 in gu minus {S} do
 eq1:=X[gu1]-add(T[gu1,gu2]*X[gu2],gu2 in gu):
 eq:=eq union {eq1}:
od:



var:=solve(eq,var):

lu:=normal(add(subs(var,X[gu1]),gu1 in gu)):

L:=normal(subs(t=1/(1-x),lu)):

factor(normal(1/(1-x/(1-x)-L))):
end:

  

#CheckGFone(C,N): checks GFone(C,x) up to N terms. Try:
#CheckGFone([2,1,2],12);
CheckGFone:=proc(C,N) local gu,i,x,mu:
gu:=GFone(C,x):
gu:=taylor(gu,x=0,N+1):
gu:=[seq(coeff(gu,x,i),i=0..N)]:
mu:=GoodCseqD(N,{C}):

if gu=mu then
 RETURN(true):
else
 RETURN(mu,gu):
fi:


end:

#Reduce11(Cset): given a set of compositions performs one step in kicking out superflous compositions. Try:
#Reduce11({[1,2,1],[3,3,3]});
Reduce11:=proc(Cset) local C1,C2:

if nops(Cset)=1 then
 RETURN(Cset):
fi:

for C1 in Cset do
 for C2 in Cset do
  if C1<>C2 and IsCon(C1,C2) then
    RETURN(Cset minus {C2}):
 fi:
 od:
od:
Cset:
end:


#Reduce1(Cset): given a set of compositions performs finds a minimal subset. Try
#Reduce1({[1,2,1],[3,3,3]});
Reduce1:=proc(Cset) local gu,gu1,gu2:
gu:=Cset:
gu1:=Reduce11(gu):

while gu<>gu1 do
gu2:=Reduce11(gu1):
gu:=gu1:
gu1:=gu2:
od:
gu:
end:




#InitialStatesS(SetC): the set of initial states for a set of compositions SetC. Try
#InitialStatesS({[1,2,1],[4,1]});
InitialStatesS:=proc(SetC) local gu,gu1:

if Reduce1(SetC)<>SetC then
 print(`Not reduced`):
 RETURN(FAIL):
fi:


gu:=powerset(SetC) minus {{}}:

{seq(Skyline(gu1),gu1 in gu)}:
end:




#InitialWt(SetC1,x,t): The  initial weight of the set SetC1. Try:
#InitialWt({[1,3,1],[3,1,3]},x,t);
InitialWt:=proc(SetC1,x,t) local gu:
option remember:
 gu:=Skyline(SetC1):

 (-1)^nops(SetC1)*t^nops(gu)*x^convert(gu,`+`):
end:

#InitialWtSt(SetC,S,x,t): The  initial weight of the state S with set of compositions SetC. Try:
#InitialWtSt({[2,1,2]},[2,1,2],x,t);
InitialWtSt:=proc(SetC,S,x,t) local lu1,gu,lu:
lu:=powerset(SetC):

if not member(S,InitialStatesS(SetC)) then
  RETURN(FAIL):
fi:

gu:=0:

for lu1 in lu do
 if Skyline(lu1)=S then
  gu:=gu+InitialWt(lu1,x,t):
fi:
od:
gu:
end:

 




#AbaSet(SetC1,S,r,x,t): inputs a set of compositions SetC1, a state S, and r between 2 and nops(Skyline(SetC)), outputs the new state and the
#change in weight in terms of x and t. Try
#AbaSet({[2,1,2,1,2]},[2,1,2,2],4,x,t);
AbaSet:=proc(SetC1,S,r,x,t) local C1,gu,mu,i:
option remember:

if Reduce1(SetC1)<>SetC1 then
 print(`Not reduced`):
 RETURN(FAIL):
fi:

C1:=Skyline(SetC1):

if ( r<=1 or r>nops(C1)) then
 print(`Bad input`):
 RETURN(FAIL):
fi:

gu:=[[op(C1),0$(r-1+nops(S)-nops(C1))],[0$(r-1),op(S)]]:
mu:=[seq(max(gu[1][i],gu[2][i]),i=1..nops(S)+r-1)]:

[[op(1..nops(S),mu)],(-1)^nops(SetC1)*t^(r-1)*x^(convert(mu,`+`)-convert(S,`+`))]:
end:




#AvotSet(SetC,S): all the fathers of the state S with a set of compositions SetC. Try:
#AvotSet({[2,1,2,1,2]},[2,1,2,2]);
AvotSet:=proc(SetC,S) local SetC1,C1,r,x,t,gu,lu:

if Reduce1(SetC)<>SetC then
 print(`Not reduced`):
 RETURN(FAIL):
fi:

gu:={}:

lu:=powerset(SetC) minus {{}}:

for SetC1 in lu do
 C1:=Skyline(SetC1):

gu:=gu union {seq(AbaSet(SetC1,S,r,x,t)[1],r=2..nops(C1))}:

od:

gu:

end:

#StatesCset(SetC): The set of states for the clusters of the set of compositions SetC. Try:
#StatesCset({[2,1,2]});
StatesCset:=proc(SetC) local lu,lu1,gu,gu1,gu2, gua:

if Reduce1(SetC)<>SetC then
 print(`Not reduced`):
 RETURN(FAIL):
fi:

lu:=InitialStatesS(SetC):
gu:=InitialStatesS(SetC):
gu1:=gu union {seq(op(AvotSet(SetC,lu1)),lu1 in lu)}:

while gu<>gu1 do
 gu2:=gu1 union {seq(op(AvotSet(SetC,gua)),gua in gu1)}:
 gu:=gu1:
 gu1:=gu2:
od:
gu1:
end:



#GFset(SetC,x): inputs a set of compositions S and outputs the generating function enumerating compositions that avoid
#the members of S as sub-compositions. Try:
#GFset({[2,3,2],[4,1,4]},x);
GFset:=proc(S,x) local SetC,gu,ku,ku1,lu,gu1,gu2,T,t,eq,var,vu,L,X,r,eq1:
option remember:
if not type(S,set) then
 print(`bad input`):
 RETURN(FAIL):
fi:

if nops({seq(nops(SetC),SetC in S)})<>1 then
 print(`The members of`, S, `are not of the same length`):
 RETURN(FAIL):
fi:


SetC:=Reduce1(S):

gu:=StatesCset(SetC):

vu:=InitialStatesS(SetC):


ku:=powerset(SetC) minus {{}}:

for  gu1 in gu do
 for gu2 in gu do
  T[gu1,gu2]:=0:
 od:
od:

for gu1 in gu do
 for ku1 in ku do
 for r from 2 to nops(Skyline(ku1)) do
  lu:=AbaSet(ku1,gu1,r,x,t):
   T[lu[1],gu1]:=T[lu[1],gu1]+lu[2]:
 od:
od:
od:



var:={seq(X[gu1],gu1 in gu)}:

eq:={}:

for gu1 in gu  do


 eq1:=X[gu1]-add(T[gu1,gu2]*X[gu2],gu2 in gu):


if member(gu1,vu) then

 eq1:=eq1-InitialWtSt(SetC,gu1,x,t):

fi:

 eq:=eq union {eq1}:
od:




var:=solve(eq,var):

lu:=normal(add(subs(var,X[gu1]),gu1 in gu)):

L:=normal(subs(t=1/(1-x),lu)):

factor(normal(1/(1-x/(1-x)-L))):
end:

#CheckGFset11(SetC,N): checks GFset(SetC,x) up to N terms. Try:
#CheckGFset11({[2,1,2]},12);
CheckGFset11:=proc(SetC,N) local gu,i,x,mu:
option remember:
gu:=GFset(SetC,x):
gu:=taylor(gu,x=0,N+1):
gu:=[seq(coeff(gu,x,i),i=0..N)]:
mu:=GoodCseqD(N,SetC):

if gu=mu then
 RETURN(true):
else
print(mu,gu):
 RETURN(false):
fi:


end:


#CheckGFset1(r,k,N): checks GFset(SetC,x) up to N terms for all non-empty subsets of the set of compositions of r with k parts. Try:
#CheckGFset1(5,2,12);
CheckGFset1:=proc(r,k,N) local SetC,gu:
gu:=powerset(Comps1(r,k)) minus {{}}:
 for SetC in gu do
  if not CheckGFset11(Reduce1(SetC),N) then
   print(`Wrong for the set`, SetC):
   RETURN(false):
 fi:
 od:
true:
end:


#CheckGFset(R,N): checks GFset(SetC,x) up to N terms for all non-empty subsets of the set of compositions of <=R with the same number of parts. Try:
#CheckGFset(3,12);
CheckGFset:=proc(R,N) local r,k:

for r from 1 to R do
 for k from 1 to r do

 if not CheckGFset1(r,k,N) then
  print(r,k):
   RETURN(false):
 fi:
od:
od:

true:

end:


#Info(S,x,N): inputs a set of compositions of the same length, SetC, and a variable x, and a positive integer N, returns
#(i) The generating function of the seqeunce enumerating compositions that avoid the compositions in SetC. let's call it a(n)
#(ii) the limit of a(n+1)/a(n), as n goes to infinity, let's call it alpha
#(iii) The limit of a(n)/alpha^n as n goes to infinity
#(iv): The first N terms of the enumerating sequence. Try:
#Info({[2,3,2]},x,30);
Info:=proc(S,x,N) local C,gu,lu,i,hal,aluf,si,alpha,gu1,Sid:
option remember:
Digits:=20:
if not type(S,set) then
 print(`bad input`):
 RETURN(FAIL):
fi:

if S={} then
 RETURN([S,normal(1+x/(1-2*x)),[1,seq(2^(i-1),i=1..N)],2,1/2]):
fi:

if nops({seq(nops(C),C in S)})<>1 then
 print(`The members of`, S, `are not of the same length`):
 RETURN(FAIL):
fi:


gu:=GFset(S,x):

gu1:=taylor(gu,x=0,max(N,201)):

Sid:=[seq(coeff(gu1,x,i),i=0..N)]:

lu:=[fsolve(denom(gu))]:

if lu=[] then
 RETURN([S,gu,Sid]):
fi:



aluf:=lu[1]:
si:=abs(lu[1]):

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

if aluf<0 then
 RETURN([S,gu,Sid]):
fi:

alpha:=1/si:



if abs(coeff(gu1,x,200)/coeff(gu1,x,199)-alpha)>10^(-12) then
 RETURN([S,gu,Sid]):
fi:


C:=-alpha/subs(x=si,diff(1/gu,x)):

if abs(coeff(gu1,x,200)/alpha^200-C)>10^(-12) then
 RETURN([S,gu,Sid,alpha]):
fi:

[S,gu,Sid,alpha,C]:
end:




#InfoV(S,x,N): Verbose form of Info(S,x,N) (q.v.). Try:
#InfoV({[2,3,2]},x,30);
InfoV:=proc(S,x,N) local gu,a,n:
gu:=Info(S,x,N):

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

print(``):
if nops(S)=1 then
print(`Theorem: Let a(n) be the number of compositions of n avoiding, as a subcomposition`, S[1]):
else
print(`Theorem: Let a(n) be the number of compositions of n avoiding as subcompositions, the members of the set`, S):
fi:
print(``):
print(`We have the following rational function for the (ordinary) generating function of a(n)`):
print(``):
print(Sum(a(n)*x^n,n=0..infinity)=gu[2]):
print(``):
print(`and in Maple format`):
print(``):
lprint(gu[2]):
print(``):
print(`The first`, N+1, `terms of a(n),  starting at n=0, are`):
print(``):
print(gu[3]):
print(``):

if nops(gu)>=4 then
print(`The limit of a(n+1)/a(n) as n goes to infinity is`):
print(``):
lprint(evalf(gu[4],12)):
print(``):
fi:

if nops(gu)=5 then
print(`a(n) is asymptotic to`):
print(``):
lprint(evalf(gu[5],12)*evalf(gu[4],12)^n):
print(``):
fi:
print(``):
print(`----------------------------------------------------------------------------`):
print(``):
end:

#InfoVth(S,x,N,co): Verbose form of Info(S,x,N) (q.v.), with theorem numbered. Try:
#InfoVth({[2,3,2]},x,30,3);
InfoVth:=proc(S,x,N,co) local gu,a,n:
gu:=Info(S,x,N):

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

print(``):
if nops(S)=1 then
print(`Theorem Number`, co, `: Let a(n) be the number of compositions of n avoiding, as a subcomposition`, S[1]):
else
print(`Theorem Number`,co, `: Let a(n) be the number of compositions of n avoiding as subcompositions, the members of the set`, S):
fi:
print(``):
print(`We have the following rational function for the (ordinary) generating function of a(n)`):
print(``):
print(Sum(a(n)*x^n,n=0..infinity)=gu[2]):
print(``):
print(`and in Maple format`):
print(``):
lprint(gu[2]):
print(``):
print(`The first`, N+1, `terms of a(n),  starting at n=0, are`):
print(``):
print(gu[3]):
print(``):

if nops(gu)>=4 then
print(`The limit of a(n+1)/a(n) as n goes to infinity is`):
print(``):
lprint(evalf(gu[4],12)):
print(``):
fi:

if nops(gu)=5 then
print(`a(n) is asymptotic to`):
print(``):
lprint(evalf(gu[5],12)*evalf(gu[4],12)^n):
print(``):
fi:
print(``):
print(`----------------------------------------------------------------------------`):
print(``):
end:


#SuperInfo1(k,x,N): Inputs a positive integer k, variable x, and positive integer N, gathers Info({C},x,N) for all
#compositions C of k. Try:
#SuperInfo1(5,x,30);
SuperInfo1:=proc(k,x,N) local lu,mu,mu1,T,lu1,ku,ru,K1,K2,i,gu:
lu:=Comps(k):
mu:={seq(GFone(lu[i],x),i=1..nops(lu))}:

for mu1 in mu do
 T[mu1]:={}:
od:

for lu1 in lu do
 ku:=GFset({lu1},x):
 T[ku]:=T[ku] union {lu1}:
od:

gu:={}:

for mu1 in mu do
 ru:=Info({T[mu1][1]},x,N):
 if nops(ru)<=3 then
 gu:=gu union {1}:
  K1[mu1]:=1:
  K2[1]:=mu1:
 else
 gu:=gu union {ru[4]}:   

   
 K1[mu1]:=ru[4]:
  K2[ru[4]]:=mu1:
   
 fi:
od:

gu:=sort(convert(gu,list)):

#[seq([gu[i],K2[gu[i]],T[K2[gu[i]]] ] ,i=1..nops(gu))]:

[seq([ gu[i],  T[K2[gu[i]]], Info({T[K2[gu[i]]][1]},x,N)] ,i=1..nops(gu))]:

end:

#SuperInfo1t(k):  Terse version of SuperInfo1(k,x,N): only gives the growth constants for all the compositions of k. Try:
#SuperInfo1t(5);
SuperInfo1t:=proc(k) local x,lu,mu,mu1,T,lu1,ku,ru,K1,K2,i,gu:
lu:=Comps(k):
mu:={seq(GFone(lu[i],x),i=1..nops(lu))}:

for mu1 in mu do
 T[mu1]:={}:
od:

for lu1 in lu do
 ku:=GFset({lu1},x):
 T[ku]:=T[ku] union {lu1}:
od:

gu:={}:

for mu1 in mu do
 ru:=Info({T[mu1][1]},x,10):
 if nops(ru)<=3 then
 gu:=gu union {1}:
  K1[mu1]:=1:
  K2[1]:=mu1:
 else
 gu:=gu union {ru[4]}:   

   
 K1[mu1]:=ru[4]:
  K2[ru[4]]:=mu1:
   
 fi:
od:

gu:=sort(convert(gu,list)):

[seq([ T[K2[gu[i]]][1],gu[i]],i=1..nops(gu))]:

end:



#SuperInfo1tV(K): Verbose form  of SuperInfo1t(k) (q.v.)  for k from 2 to K: Try:
#SuperInfo1tV(5);
SuperInfo1tV:=proc(K) local gu,i,k,t0:

t0:=time():

print(`All the growth constants of Compositions that add up to at most`, K):
print(``):
print(`By Shalsoh B. Ekhad `):
print(``):
print(`In this article we will list all the growth constants, in ascending order, of the enumnerating sequences`):
print(`avoding, by containment, the compsoitions that sum-up to less than`, K):
print(``):
if K<=9 then
print(`for readabiltiy, we will get rid of the commas`):
else
print(`for readabiltiy, we will get rid of the commas when the parts are <=9`):
fi:

for k from 2 to min(K,9) do
gu:=SuperInfo1t(k):
 print(`The growth constants, for compositions of `, k, `in ascending order (where we only show one representative for equivalent compositions) are`):

for i from 1 to nops(gu) do
print(cat(op(gu[i][1])) ,  gu[i][2]):
od:
print(`-----------------------------------------------`):
od:

print(``):
for k from 10 to K  do
gu:=SuperInfo1t(k):
 print(`The growth constants, for compositions of `, k, `in ascending order (where we only show one represenative for equivalent compositions) are`):

for i from 1 to nops(gu) do
print(op(gu[i][1]), `:`, gu[i][2]):
od:

print(`-----------------------------------------------`):

od:

print(``):
print(`-----------------------------------------------------------`):
print(``):
print(`This ends this article, that took`, time()-t0, `to create. `):

end:



#SuperInfo1V(k,x,N): verbose verson of SuperInfo1(k,x,N) (q.v.). Try:
#SuperInfo1V(5,x,30);
SuperInfo1V:=proc(k,x,N) local gu,co,t0:
t0:=time():
gu:=SuperInfo1(k,x,N):



print(``):
print(`----------------------------------------`):
print(``):

print(`Generating functions and Growth rates for the Enumerating Sequences for Comositions Avoiding Each of the compositions of`, k):
print(``):
print(`By Shalosh B. Ekhad`):
print(``):
co:=1:
print(`The compositions of`, k, `that yield polynomial growth (hence growth ratio 1, that is the lowest) are`, op(gu[1][2])):
print(``):
InfoVth({gu[1][2][1]},x,N,co):
print(``):

for co from 2 to nops(gu) do
print(`The compositions of`, k, `that yield the`,co, `-th largest growth, that is`, gu[co][1], ` are `,  op(gu[co][2]) ):
print(``):
InfoVth({gu[co][2][1]},x,N,co):
print(``):

od:
print(``):
print(`This ends this article, that took`, time()-t0, `seconds to generate. `):
print(``):
print(`----------------------------------------------------------`):
print(``):
end:



#SuperInfo2kr(k,r,x,N): Inputs  positive integers k and r, variable x, and positive integer N, gathers Info(S,x,N) 
#for all subsets of the set of compositions of k with r parts. Try:
#SuperInfo2kr(5,2,x,30);
SuperInfo2kr:=proc(k,r,x,N) local lu,mu,mu1,T,lu1,ku,ru,K1,K2,i,gu:

lu:=Comps1(k,r):
lu:=powerset(lu) minus {{}}:

mu:={seq(GFset(lu[i],x),i=1..nops(lu))}:


for mu1 in mu do
 T[mu1]:={}:
od:

for lu1 in lu do
 ku:=GFset(lu1,x):
 T[ku]:=T[ku] union {lu1}:
od:


gu:={}:

for mu1 in mu do
 ru:=Info(T[mu1][1],x,N):
 if nops(ru)<=3 then
 gu:=gu union {1}:
  K1[mu1]:=1:
  K2[1]:=mu1:
 else
 gu:=gu union {ru[4]}:   

   
 K1[mu1]:=ru[4]:
  K2[ru[4]]:=mu1:
   
 fi:
od:

gu:=sort(convert(gu,list)):



#[seq([gu[i],K2[gu[i]],T[K2[gu[i]]] ] ,i=1..nops(gu))]:

[seq([ gu[i],  T[K2[gu[i]]], Info(T[K2[gu[i]]][1],x,N)] ,i=1..nops(gu))]:

end:


#SuperInfo2krV(k,r,x,N): verbose verson of SuperInfo2kr(k,r,x,N) (q.v.). Try:
#SuperInfo2krV(5,2,x,30);
SuperInfo2krV:=proc(k,r,x,N) local gu,co,t0:
t0:=time():
gu:=SuperInfo2kr(k,r,x,N):



print(``):
print(`----------------------------------------`):
print(``):

print(`Generating functions and Growth rates for the Enumerating Sequences for Comositions Avoiding Each of the Possible subsets of`):
print(`the set of compositions of`, k, `into exactly`, r, `parts `):
print(``):
print(`By Shalosh B. Ekhad`):
print(``):
for co from 1 to nops(gu) do
print(`The subsets of the set of compositions of`, k, `with exactly`, r, `parts, that yield the`):
print(co, `-th largest growth, that is`, gu[co][1], ` are `,  op(gu[co][2]) ):
print(``):
InfoVth(gu[co][2][1],x,N,co):
print(``):

od:
print(``):
print(`This ends this article, that took`, time()-t0, `seconds to generate. `):
print(``):
print(`----------------------------------------------------------`):
print(``):
end: