#SmallExp.txt: One of the Maple programs accompanying the article
#"Counting Permutations that Avoid Many Patterns" by
#Yonah BIERS-ARIEL, Haripriya CHAKRABORTY, John CHIARELLI, Bryan EK, Andrew LOHR, 
#Jinyoung PARK, Justin SEMONSEN, Richard VOEPEL, Mingjia YANG, Anthony ZALESKI, 
#and Doron ZEILBERGER 
#Available from:
#http://www.math.rutgers.edu/~zeilberg/mamarim/mamarimhtml/pamp.html
#JohnChiarelliFullCode.txt.txt: One of the Maple programs accompanying the article
#"Counting Permutations that Avoid Many Patterns" by
#Yonah BIERS-ARIEL, Haripriya CHAKRABORTY, John CHIARELLI, Bryan EK, Andrew LOHR, 
#Jinyoung PARK, Justin SEMONSEN, Richard VOEPEL, Mingjia YANG, Anthony ZALESKI, 
#and Doron ZEILBERGER 
#Available from:
#http://www.math.rutgers.edu/~zeilberg/mamarim/mamarimhtml/pamp.html

# Ask(num,N,n): generates a random set S of num permutations of
# length 4, finds sequence of number of permutations of length n
# avoiding S for n from 1 to N; outputs low-degree polynomial that
# sequence matches and point at which it starts doing so if it
# exists, otherwise outputs FAIL and sequence itself
Ask:=proc(num,N,n) local S,per,te,sol,i,k,oh:
 S:={}:
 per:=[op(GG(4,{}))]:
 while nops(S)<num do
  te:=rand(1..24)():
  S:=S union {per[te]}:
 od:
# print(S):
 sol:=SeqGG(S,N):
# print(sol):
 for k from 1 to N-2 do
  oh:=GP([[seq(i,i=k..N)],[seq(sol[i],i=k..N)]],n):
  if oh<>FAIL then
   RETURN([oh,k]):
  fi:
 od:
 [FAIL,sol]:
end:  

# Receive(num,N,T,n): runs Ask(num,N,n) T times, counts number of
# times resulting sequence goes to 0 or fails to match a polynomial,
# and how many of the rest have each possible degree
Receive:=proc(num,N,T,n) local zed,pol,big,i,quer:
 zed:=0:
 pol:=[0$(N-1)]:
 big:=0:
 for i from 1 to T do
  quer:=Ask(num,N,n):
  if quer[1]=FAIL then
   big:=big+1:
   print(quer[2]):
  elif quer[1]=0 then
   zed:=zed+1:
  else
   quer:=degree(quer[1],n)+1:
   pol[quer]:=pol[quer]+1:
  fi:
 od:
 print(`We tested `,T,` sets of `,num,` permutations for the number of satisfying sequences of length n.`):
 print(zed,` of them went to 0 in the limit.  Of the remainder:`):
 for i from 1 to N-1 do
  if pol[i]>0 then
   print(pol[i],` of them approached a polynomial of degree `,i-1):
  fi:
 od:
 print(big,` of them did not approach a polynomial within `,N,` steps.`):
end:

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

#C27.txt,  Dec. 8, 2016 Pattern Avoiding Permutations
Help:=proc(): print(` redu(L), IV(n,k) , IsBad(pi,P) , GG(n,P) , SeqGG(P,N), GP1(L.n,d), GP(L,n) `): end:

with(combinat):

#inputs a list of numbers, and outputs the reduction as a permutation of nops(L)
#redu([e,pi,gamma,phi, h])=[ 4, 5 , 2 , 3    ,1]
redu:=proc(L) local L1,T,k,i:
L1:=sort(sort(L,'output'='permutation'),'output'='permutation'):

end:

#IV(n,k): the set of increasing vectors of length k in {1, ...,n} ending in n
IV:=proc(n,k) local S,i,s:
S:=choose([seq(i,i=1..n-1)],k-1):
{seq([op(s), n], s in S)}:
end:

 #IsBad1(pi,p): inputs a perm. pi and a pattern p outputs true iff pi contains the pattern p
 #where the last entry of pi participates
 IsBad1:=proc(pi,p) local S,n,k,s:
  n:=nops(pi):
  k:=nops(p):
  S:=IV(n,k):
  for s in S do
   if redu(pi[s])=p then
      RETURN(true):
   fi:
  od:

 false: 
   end:
 #inputs a permutation pi and a set of patterns P and outputs true if
 #it contains at least one member of P
 IsBad:=proc(pi,P) local p:
 for p in P do
  if IsBad1(pi,p) then
     RETURN(true):
   fi:
 od:
  false:
 end:
   
   #inputs a non-neg. integer n and a set of patterns P and outputs the
   #set of permutations avoiding the patterns of P 
   GG:=proc(n,P) local Sold,S, i,cand,pi:
   option remember:
    if n=0 then
       RETURN({[]}):
    fi:
    Sold:=GG(n-1,P):
    S:={}:

    for pi in Sold do
     for i from 0 to n-1 do
     cand:=redu([op(pi),i+1/2]):
     if not IsBad(cand,P) then
        S:=S union {cand}:
     fi:
     od:
   od:
    S:
     end:
  
   #SeqGG(P,N): inputs a set of patterns P, and a positive integer N
    #returns the list of length N whose i-th entry is the NUMBER of permutations
    #of length i NOT containing any of the patterns in P
   SeqGG:=proc(P,N) local i:
    [seq(nops(GG(i,P)),i=1..N)]:
   end:








#GP(L,n): inputs a list of pairs [input,output], and a variable n
#and outputs a polynomial of degree <=nops(L) such P(input)=output
GP:=proc(L,n) local d,hope:
for d from 0 to nops(L[1])-3 do
 hope:=GP1(L,n,d):
# print(d,hope):
 if hope<>FAIL then
   RETURN(hope):
 fi:
od:
FAIL:
end:

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

 
#GP1(L,n,d): inputs a list of pairs [input,output], and a variable n
#and outputs a polynomial of degree d such P(input)=output
GP1:=proc(L,n,d) local lo,P,i,a,var,eq:
#print(op(L)):
lo:=L[1]:
if nops(lo)<d+2 then
 RETURN(FAIL):
fi:
P:=add(a[i]*n^i,i=0..d):
#print(P):
var:={seq(a[i],i=0..d)}:
eq:={seq(subs(n=L[1][i],P)=L[2][i],i=1..d+2)}:
#print(eq):
var:=solve(eq,var):
#print(var):
if var=NULL then
  RETURN(FAIL):
fi:

P:=subs(var,P):
#print(P):
#print(seq(subs(n=L[1][i],P)-L[2][i],i=1..nops(lo))):
if {seq(subs(n=L[1][i],P)-L[2][i],i=1..nops(lo))}<>{0} then
# print(fie):
  RETURN(FAIL):
else
 RETURN(P):
fi:

end: