######################################################################
##RPneg.txt: Save this file as  RPneg.txt                            #
## To use it, stay in the                                            #
##same directory, get into Maple (by typing: maple <Enter> )         #
##and then type:  read RPneg.txt<Enter>                              #
##Then follow the instructions given there                           #
##                                                                   #
##Written by Mingjia Yang, Rutgers University                        #
#my237 at math dot rutgers dot edu                                   #
######################################################################
 
#Created: May 2019
 
print(`Created:  May 2019`):
print(` This is RPneg.txt `):
print(`It is one of the packages that accompany the article `):
print(` Systematic Counting of Restricted Partitions`):
print(`by Mingjia Yang and Doron Zeilberger`):

print(`and also 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/  .`):

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



Help:=proc():

if args=NULL then
 print(`The main procedures are: checkP,C,P,PKM,Pkm, PB , xnSeqNeg `):
 print(` `):


elif nops([args])=1 and op(1,[args])=PKM then
print(`PKM(A,K,M,q):input integers K and M, and a set of patterns A, output the generating function for partitions starting with a number at most K with at most M parts and avoid patterns in A, using the cluster method`):
print(`Try: `):
print(`PKM({[0],[2,5]},10,10,q); `):


elif nops([args])=1 and op(1,[args])=Pkm then
print(`Pkm(A,k,m,q): input integers k and m, and a set of patterns A, output the generating function for partitions starting with k with length m and avoid patterns in A, using the cluster method`):
print(`Try: `):
print(`Pkm({[0],[2,5]},10,5,q); `):


elif nops([args])=1 and op(1,[args])=C then
print(`C(A,k,l,w,q):the weight of clusters starting with k, ending with l, with width w, with the set of forbidden patterns A`):
print(`Try: `):
print(`C({[0],[1,1]},4,1,4,q); `):

elif nops([args])=1 and op(1,[args])=P then
print(`P(A,m,q):input an integer m and a set of patterns A, output the first m terms of the generating function for partitions that avoid patterns in A, using the cluster method`):
print(`Try: `):
print(`P({[1,1],[0]},10,q); `):


elif nops([args])=1 and op(1,[args])=checkP then
print(`checkP(A,m,q):checks P up to the mth term with the set of patterns A`):
print(`Try: `):
print(`checkP({[0],[1,3,5],[5,6,7,8]},10,q); `):

elif nops([args])=1 and op(1,[args])=PB then
print(`PB(m,A,q):input an integer m and a set of patterns A, output the first m terms of the generating function for partitions that avoid patterns in A`):
print(`Try: `):
print(`PB(10,{[1,1],[0]},q);`):


elif nops([args])=1 and op(1,[args])=xnSeqNeg then
print(`xnSeqNeg(m,A): The first  m terms of the sequence enumerating partitions that avoid the set of patterns in A. Try:`):
print(`xnSeqNeg(12,{[0,0],[0,1]});`):


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

#################################################

overlap:=proc(u,v,q) local i,j,k,lu,gug:
 option remember:
lu:=[]:
 
for i from 2 to nops(u) do
  for j from i to nops(u) while (j-i+1<=nops(v) and op(j,u)=op(j-i+1,v))
         do
   :
  od:   
 
 if j=nops(u)+1 and (i>1 or j>2) then
  gug:=1:
 
  for k from 1 to i-1 do
   gug:=gug*q^op(k,u):
  od:
 
  lu:=[op(lu),[gug,nops(u)-i+1]]:
 
fi:
 
od:
 
lu:
 
end:


overlap1:=proc(u,v,k1,k2,q) local u1, v1,i,j:
option remember:
u1:=[k1,seq(k1-add(u[i],i=1..j),j=1..nops(u))]:
v1:=[k2,seq(k2-add(v[i],i=1..j),j=1..nops(v))]:

overlap(u1,v1,q):
end:



#issubword(u,v) returns 1 if v is a subword of u, otherwise 0
 
issubword:=proc(u,v) local i,j:
option remember:

 
 
for i from 1 to nops(u) do
  for j from i to nops(u) while (j-i+1<=nops(v) and op(j,u)=op(j-i+1,v))
         do
   :
  od:   
 
if j-i=nops(v)  then
 RETURN(1):
 fi:
od:
 
0:
 
end:

#superflous(B,w) given a set of bad words, and one of its
#members, finds whether it is superflous,
#that is, whether it is a factor of another one, and deletes the larger one
 
superflous:=proc(B,w)
local i,v:
 
if not type(B,set) or not type(w,list) then
 ERROR(`Bad input`):
fi:
 
if not member(w,B) then
 ERROR(`Second argument must belong to first arg.`):
fi:
 
for i from 1 to nops(B) do
v:=op(i,B):
if v<>w and issubword(w,v)=1 then
 RETURN(1):
fi:
od:
0:
end:

hakten:=proc(B)
local w,i:
 
for i from 1 to nops(B) do
w:=op(i,B):
 
if superflous(B,w)=1 then
 RETURN(B minus {w}):
fi:
od:
B:
end:
 
#Hakten(B) removes all superflous words
 
Hakten:=proc(B)
local B1,B2:
 
B1:=B:
 
B2:=hakten(B1):
 
while B1<>B2 do
B1:=B2:
B2:=hakten(B2):
od:
 
B2:
end:

####################################

C1:=proc(A,v,k,l,w,q) local g,i,k1,p,u:
option remember:

if k<l then
RETURN(0):
fi:

if w<0 then
RETURN(0):
fi:

g:=-q^(k+add(k-add(v[i],i=1..j),j=1..nops(v)))*`if`(k=l+add(v[i],i=1..nops(v)) and w=nops(v)+1,1,0):


for k1 from 1 to k do
for u in A do
for p in overlap1(v,u,k,k1,q) do
g:=g-p[1]*C1(A,u,k1,l,w-nops(v)-1+p[2],q):
od:
od:
od:

g:

end:

C:=proc(A,k,l,w,q) option remember:

add(C1(A,v,k,l,w,q),v in A):
end:

######
Pkm:=proc(A,k,m,q) local l,r,w,i,a,A1:
option remember:

A1:=Hakten(A):

if k=1 then
for a in A1 do 
if convert(a,set)={0} and m>nops(a) then
RETURN(0):
fi:
od:
RETURN(q^m):
fi:


if m<0 then
RETURN(0):
fi:

if m=0 then
RETURN(1):
fi:

if m=1 then
RETURN(q^k):
else
expand(q^k*add(Pkm(A1,r,m-1,q),r=1..k)+add(add(`if`(m=w,1,add(Pkm(A1,r,m-w,q),r=1..l))*C(A1,k,l,w,q),l=1..k),w=1..m)):

fi:
end:


P:=proc(A,m,q) option remember:
convert(taylor(PKM(A,m,m,q),q=0,m+1),polynom):
end:

checkP:=proc(A,m,q):
if P(A,m,q)-PB(m,A,q)=0 then
RETURN(true):
fi:
false:
end:


PKM:=proc(A,K,M,q) option remember:

add(add(Pkm(A,k,m,q),m=1..M),k=1..K):

end:


 

############################################bruce force below#####
#input an integer m, output the first m terms of the generating function for partitions that avoid patterns in A
PB:=proc(m,A,q) local n:

add(NP(n,A)*q^n,n=1..m):

end:

#coefficients of PB

CPB:=proc(m,A,q) local i:

seq(coeff(PB(m,A,q),q,i),i=0..m):

end:

#coefficients of P

CP:=proc(m,q) local i:

seq(ifactor(coeff(P(m,q),q,i)),i=0..m):

end:





#NP(n,A) inputs an integer n and a set of patterns A, output the number of partitions of n
#that avoid patterns in A.

NP:=proc(n,A):
nops(ParA(n,A)):
end:


#set of partitions of n with the largest part k

Par:=proc(n,k) local p,k1:
option remember:
if n=k then
RETURN({[k]}):
fi:
if n<k then
RETURN({}):
fi:
`union`(seq(seq({[k,op(p)]},p in Par(n-k,k1)),k1=1..k)):
end:


#ParA takes a set of patterns A and an integer n and output the set of 
#partitions of n that does not contain patterns in A

ParA:=proc(n,A) local k: 
`union`(seq(PA(A,n,k),k=1..n)):
end:



#PA takes a set of patterns A and integers k and n and output the set of 
#partitions of n with largest part k and does not contain patterns in A

PA:=proc(A,n,k) local p,S:
S:={}:
for p in Par(n,k) do 
if not Ct(p,A) then
S:=S union {p}:
fi:
od:
S:
end:


#Ct(p,A) takes a partition p and examines whether it contains the patterns in A 

Ct:=proc(p,A) local a:
for a in A do
if Ct1(p,a) then
RETURN(true):
fi:
od:
false:
end:


#Ct1(p,a) takes a partition p and examines if it "contains" the pattern a (for example, if a=[1,2], 
#then if p[i]-p[i-1]=1 and p[i-1]-p[i-2]=2 for some i, then p "contains" the pattern a)

Ct1:=proc(p,a) local m,r,i:
m:=nops(a):
for i from 1 to nops(p)-m do
if {seq(evalb(p[i+r-1]-p[i+r]=a[r]),r=1..m)}={true} then
RETURN(true): 
fi:
od:
false:
end:





###########################Another version of the brute force, using reccursion#########################

#ParA1(M,N,A): inputs two integers M and N, output the set of partitions with the first part<=M and no more 
#than N parts and avoid patterns in A 

ParA1:=proc(M,N,A) local S1, S2, p,M1:
option remember:

if N=1 then
RETURN({seq([M1],M1=1..M)})
fi:

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


S2:={}:

S1:=`union`(seq(ParA1(M1,N-1,A),M1=1..M)): 

for p in S1 do
if not Ct([M,op(p)],A) then
S2:=S2 union {[M,op(p)]}:
fi:
od:

S2 union ParA1(M-1,N,A):

end:

ParA2:=proc(M,N,A) local N1:
option remember:

`union`(seq(ParA1(M,N1,A), N1=1..N)):

end:


#Input an interger M, a set of patterns A, output the first M terms of the generating function for partitions avoiding 
#patterns in A

ParA2:=proc(M,A,q) local p,g,i:

option remember:

g:=0:

for p in ParA2(M,M,A) do
g:=g+q^add(p):
od:

add(coeff(g,q^i)*q^i,i=1..M):


end:



CPB1:=proc(m,A,q) local i:
option remember:

seq(coeff(ParA2(m,A,q),q,i),i=0..m):

end:




#xnSeqNeg(m,A): The first  m terms of the sequence enumerating partitions that avoid the set of patterns in A. Try:
#xnSeq(20,{[0,0],[0,1]});
xnSeqNeg:=proc(m,A) local gu,q,i:
gu:=P(A,m,q):
[seq(coeff(gu,q,i),i=1..m)]:

end: