print(` Version of April 18, 1997`):
print(`This is SYMGJ, A Maple package`):
print(`Written by John Noonan and Doron Zeilberger.`):
lprint(``):
print(`It is one of the packages accompanying Noonan and Zeilberger's `):
print(`article: "The Goulden-Jackson Cluster Method: Extensions,`):
print(`Applications, and Implementations", where the method`):
print(`is explained in detail, and then extended and generalized in`):
print(`various directions. `):
lprint(``):
print(`The paper, and the packages, including the present one`):
print(`is available from http://www.math.temple.edu/~zeilberg/gj.html`):
lprint(``):
print(`This package implements the case where the set of mistakes`):
print(`is symmetric with respect to the symmetric group on the alphabet`):
lprint(``):
print(`For general help, and a list of the available functions,`):
print(` type "ezra();". For specific help type "ezra(procedure_name)" `):
lprint(``):
 
lprint(`Warning: t,s,C,x are global variables, do not use them!`):
 
 
ezra:=proc()
if args=NULL then
print(`SYMGJ`):
print(`A Maple package that implements the Goulden-Jacskon Cluster method`):
print(` for finding generating functions and series expansions`):
print(`for the number of words avoiding a prescribed, finite, set of`):
print(`"mistakes". Contains Procedures `):
print(`SGJseries_slow,SGJst,SGJs,SGJseries,Ssolser `):
print(` For specific help type "ezra(procedure_name)" `):
fi:
 
 
if nops([args])=1 and op(1,[args])=`Ssolser` then
print(`Ssolser(vars,pols,equs,NUM,L): Given a system of eqations of the form`):
print(`a[i](s)=p[i](s)+ \sum_{j=1}^{n} P[i,j](s) a[j](s), finds the series `):
print(`expansion of sum_{i \in var} (NUM-zm)!/zm!*a[i].`):
print(`The inputs are the list of a[i] , vars, the corresponding`):
print(` list pols, of the p[i]'s, the list of expressions`):
print(` \sum_{j=1}^{n} P[i,j](s) a[j](s),  MUM is the number of letters `):
print(`subset T, and the expansion should be up to the L's term`):
fi:
 
 
if nops([args])=1 and op(1,[args])=`SGJseries` then
print(`SGJseries(NUMLETTERS,MISTAKES,L) finds the first L terms `):
print(`of the sequence: Number of words of length n, having no factors`):
print(`(i.e. subsequences of consecutive letters) that belong to the`):
print(`list of lists MISTAKES, where the number of letters of the alphabet`):
print(` is NUMLETTERS. For example to find the sequence a_n (0<=n<=15)`):
print(`where a_n:= number of words in {1,2,3} avoiding ijk (i<>j<>k)`):
print(` type SGJseries(3,[[1,2,3]],15)`):
fi:
 
 
if nops([args])=1 and op(1,[args])=`SGJst` then
print(`SGJst(NUMLETTERS,MISTAKES) finds the generating function`):
print(`F(s;t), the weight enumerator  `):
print(`of all words in the alphabet alphabet=[1,NUMLETTERS], where the`):
print(`the weight of a word is the product of s^length times`):
print(`t^(#number of factors that belong to the set MISTAKES`):
print(`and its images under the symmetric group on {1, ...NUMLETTERS}`):
print(`the set of mistakes, MISTAKES should be such that if it contains`):
print(` a word, it cannot contain any of its proper subwords`):
print(`For example to get the g.f. for the set of words`):
print(`in the alphabet {1,2,3}, avoiding 123,132,213,231,312,231, type`):
print(` type SGJst(3,{[1,2,3]});`):
fi:
 
 
 
 
if nops([args])=1 and op(1,[args])=`SGJs` then
print(`SGJs(NUMLETTERS,MISTAKES) finds the generating function`):
print(`g(s), is the g.f. sum(a[i]*s^i,i=0..infty)  `):
print(`where a[i]:=number of words of length i without any mistakes`):
print(`from MISTAKES and its images`):
print(`For example to get the g.f. for a_n:=number of words of length n`):
print(`in the alphabet {1,2,3}, avoiding 123,132,231,213,312, and 321 type`):
print(` type SGJs(3,[[1,2,3]])`):
fi:
 
end:
 
 
#overlapz is a procedure that given two words u and v
#computes the weight-enumerator of all v\suffix(u),
#for all suffixes of u that are prefixes of v, but with uniform weight s
overlapz:=proc(u,v)
local i,j,k,lu,gug:
 
lu:=0:
 
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-i=nops(v) and u<>v then
 ERROR(v,`is a subword of`,u,`illegal input`):
 fi:
 
 
 
 if j=nops(u)+1 and (i>1 or j>2) then
  gug:=1:
 
  for k from j-i+1 to nops(v) do
   gug:=gug*s:
  od:
 
  lu:=lu+gug:
 
fi:
 
od:
 
lu:
 
end:
 
 
 
 
#findeqz sets up the equ C[v]= s+t*Sum_u overlap(u,v) *C[u]
 
findeqz:=proc(v,MISTAKES)
local eq,i,u:
eq:=t-1:
 
for i from 1 to nops(v) do
 eq:=eq*s:
od:
 
for i from 1 to nops(MISTAKES) do
 u:=op(i,MISTAKES):
 eq:=eq+(t-1)*overlapz(u,v)*C[op(canon(u))]:
 
od:
 
C[op(v)]-eq=0:
 
end:
 
 
 
#findeqz_fast sets up the equ C[v]= s+t*Sum_u overlap(u,v) *C[u]
#just like findeqz but returns seperately the polynomial part
#and the equation part
 
findeqz_fast:=proc(v,MISTAKES)
local eq,PO,i,u:
 
 PO:=(-1)*s^nops(v):
 
eq:=0: 
for i from 1 to nops(MISTAKES) do
 u:=op(i,MISTAKES):
 eq:=eq-overlapz(u,v)*C[op(canon(u))]:
od:
 
PO,eq:
 
end:
 
 
SGJst:=proc(NUMLETTERS,MISTAKES1)
local v,eq, var,i,lu,MISTAKES,zm,MISTAKESFULL:
 
MISTAKES:={}:
 
for i from 1 to nops(MISTAKES1) do
MISTAKES:=MISTAKES union {canon(op(i,MISTAKES1))}:
od:
 
MISTAKESFULL:={}:
 
for i from 1 to nops(MISTAKES) do
MISTAKESFULL:=MISTAKESFULL union khaverim(op(i,MISTAKES),NUMLETTERS):
od:
 
eq:={}:
var:={}:
 
for i from 1 to nops(MISTAKES) do
 v:=op(i,MISTAKES):
 eq:= eq union {findeqz(v,MISTAKESFULL)}:
 var:=var union {C[op(v)]}:
od:
 
var:=solve(eq,var):
 
lu:=1:
 
 lu:=lu-s*NUMLETTERS:
 
for i from 1 to nops(MISTAKES) do
 v:=op(i,MISTAKES):
  zm:=nops(convert(v,set)):
 lu:=lu-NUMLETTERS!/(NUMLETTERS-zm)!*subs(var,C[op(v)]):
od:
 
normal(1/lu):
end:
 
 
 
 
SGJs:=proc(NUMLETTERS,MISTAKES1)
local v,eq, var,i,lu,MISTAKES,zm,MISTAKESFULL:
 
MISTAKES:={}:
 
for i from 1 to nops(MISTAKES1) do
MISTAKES:=MISTAKES union {canon(op(i,MISTAKES1))}:
od:
 
MISTAKESFULL:={}:
 
for i from 1 to nops(MISTAKES) do
MISTAKESFULL:=MISTAKESFULL union khaverim(op(i,MISTAKES),NUMLETTERS):
od:
 
eq:={}:
var:={}:
 
for i from 1 to nops(MISTAKES) do
 v:=op(i,MISTAKES):
 eq:= eq union {findeqz(v,MISTAKESFULL)}:
 var:=var union {C[op(v)]}:
od:
 
eq:=subs(t=0,eq):
var:=solve(eq,var):
 
lu:=1:
 
 lu:=lu-s*NUMLETTERS:
 
for i from 1 to nops(MISTAKES) do
 v:=op(i,MISTAKES):
  zm:=nops(convert(v,set)):
 lu:=lu-NUMLETTERS!/(NUMLETTERS-zm)!*subs(var,C[op(v)]):
od:
 
normal(1/lu):
end:
 
 
 
 
SGJseries:=proc(NUMLETTERS,MISTAKES1,L)
local v,eq, var,i,lu,mu,pols,ku,MISTAKES,MISTAKESFULL:
eq:=[]:
var:=[]:
pols:=[]:
 
MISTAKES:={}:
 
for i from 1 to nops(MISTAKES1) do
MISTAKES:=MISTAKES union {canon(op(i,MISTAKES1))}:
od:
 
 
 
MISTAKESFULL:={}:
 
for i from 1 to nops(MISTAKES) do
MISTAKESFULL:=MISTAKESFULL union khaverim(op(i,MISTAKES),NUMLETTERS):
od:
 
 
for i from 1 to nops(MISTAKES) do
 v:=op(i,MISTAKES):
 mu:=findeqz_fast(v,MISTAKESFULL):
 eq:= [op(eq),mu[2]]:
 pols:= [op(pols),mu[1]]:
 var:=[op(var),C[op(v)]]:
od:
 
 
lu:=Ssolser(var,pols,eq,NUMLETTERS,L+2):
 
 
ku:=1-NUMLETTERS*s:
 
for i from 1 to L+2 do
 ku:=ku-op(i,lu)*s^i:
od:
 
ku:=taylor(1/ku,s=0,L+1):
 
lu:=[]:
 
for i from 0 to L do
lu:=[op(lu),coeff(ku,s,i)]:
od:
 
lu:
end:
 
 
 
 
 
#given a word, computes all its shifts
 
shift1:=proc(u)
local gu,i:
gu:=[u]:
 
for i from 2 to nops(u) do
 gu:=[op(gu),[op(i..nops(u),u),op(1..i-1,u)]]:
od:
 
gu:
 
end:
 
 
 
 
#overlap1 finds the list of places i in v s.t. v[1],...,v[i] is a suffix of u
overlap1:=proc(u,v)
local i,gu:
gu:=[]:
 
for i from 1 to nops(v)-1 do
 if nops(u)-i+1>1 and [op(1..i,v)]=[op(nops(u)-i+1..nops(u),u)] then
  gu:=[op(gu),nops(v)-i]:
 fi:
od:
 
gu:
 
end:
 
 
 
#Given a set of mistakes overlapmat1 is a list of lists A such that 
#op(i,op(j,A)) is the overlap(op(j,MISTAKES),op(i,MISTAKES)))
 
overlapmat1:=proc(MISTAKES)
local gu,u,v,lu,i,j:
gu:=[]:
 
for i from 1 to nops(MISTAKES) do
v:=op(i,MISTAKES):
lu:=[]:
 
  for j from 1 to nops(MISTAKES) do
   u:=op(j,MISTAKES):
    lu:=[op(lu),overlap1(u,v)]:
  od:
 
 gu:=[op(gu),lu]:
od:
 
gu:
 
end:
 
 
 
Ssolser:=proc(vars,pols,equs,NUM,L)
local Mem,i,aj,maxi,X,r,lu,zm,gu, j, lu1, n, pij, v:
 
#Mem is a list of the memories that each a[j] needs
 
 
Mem:=[]:
for j from 1 to nops(vars) do
aj:=op(j,vars):
 
maxi:=0:
 
for i from 1 to nops(equs) do
 if degree(coeff(op(i,equs),aj,1),s)>maxi then
     maxi:=degree(coeff(op(i,equs),aj,1),s):
 fi:
od:
 
Mem:=[op(Mem),maxi]:
od:
 
 
for i from 1 to nops(vars) do
X[i]:=[seq(0,j=1..op(i,Mem))]:
od:
 
lu:=[]:
 
for n from 1 to L do
 for i from 1 to nops(equs) do
   gu[i]:=coeff(op(i,pols),s,n):
 
  for j from 1 to nops(vars) do
      aj:=op(j,vars):
      pij:=coeff(op(i,equs),aj,1):
 
       for r from 1 to degree(pij,s) do
         gu[i]:=gu[i]+coeff(pij,s,r)*op(nops(X[j])-r+1,X[j]):
       od:
  od:
 od:
 
lu1:=0:
 
 
for i from 1 to nops(vars) do
 v:=[op(op(i,vars))]:
 zm:=nops(convert(v,set)):
    lu1:=lu1+NUM!/(NUM-zm)!*gu[i]:
od:
 
lu:=[op(lu),lu1]:
 
for i from 1 to nops(vars) do
 X[i]:=[op(2..nops(X[i]),X[i]),gu[i]]:
od:
 
od:
lu:       
   
end:
 
 
with(combinat):
 
#canon(mila):Given a word in an alphabet of integers, in terms of a list
#finds an image under the action of the symmetric group that
#is least lexicographically
 
canon:=proc(mila)
local i,gu,mu,lu:
 
gu:=[]:
 
for i from 1 to nops(mila) do
 if nops({op(1..i,mila)})>nops({op(1..i-1,mila)}) then
    gu:=[op(gu),op(i,mila)]:
 fi:
od:
 
for i from 1 to nops(gu) do
 mu[op(i,gu)]:=i:
od:
 
lu:=[]:
 
for i from 1 to nops(mila) do
lu:=[op(lu),mu[op(i,mila)]]:
od:
 
lu:
end:
 
#peula(mila,perm): Given a word in the alphabet {1,2,..,r} applies
#to it the (partial) permutation perm i.e. i is replaced by perm[i]
 
peula:=proc(mila,perm)
local gu,i:
gu:=[]:
 
for i from 1 to nops(mila) do
gu:=[op(gu),op(op(i,mila),perm)]:
od:
 
gu:
end:
 
 
#khaverim(mila,n):Given a word mila finds the set of all its n(n-1)..(n-j+1)
#images under the symmetric group on {1,2, ..n},where j is the number of
#distinct letters in mila
 
khaverim:=proc(mila1,n)
local i,gu,r,lu,mila:
r:=nops(convert(mila1,set)):
mila:=canon(mila1):
 
if n<r then
 ERROR(`Not enough letters`):
fi:
 
gu:=permute(n,r):
 
lu:={}:
 
for i from 1 to nops(gu) do
lu:=lu union {peula(mila,op(i,gu))}:
od:
 
lu:
 
end: