# Homework 19.  April-2-2012.  Pat Devlin #
# I have no preferrence about keeping my work private. #


Help:=proc():
  print(`SFP(d, n, blank:=B), NdnC(d,n)`):
  print(`WnX(A, FP, B, n), WnGX(A,FP,B,n,x)`):
  print(`adnC(d,n)`):
  print(`WnR(A, FP, B, N)`):
  print(`Type HelpAll() for all functions`):
end proc:


HelpAll:=proc():
  print(`Adn(d,n), IsBad(d,S),adnS(d,n)`):
  print(`NdnS(d,n) , AdnM(d,n), adnM(d,n), NdnM(d,n) `):
  print(`GuessRec(L)`):
  print(`AdjustBFP(a,B,FP,BFP) , Wn(A,FP,B,n), WGn(A,FP,B,n,BFP)`):
end:


#############
# Problem 1 #
#############

# The forbidden patterns are exactly those forbidden sets found in homework 18::problem 4, reproduced below
# The choice of FP to make it so that Wn({0,1}, FP, B, n) is the same as
# Ndn(d,n)  is just avoiding all patterns of the form
# [1, B$k, 1, B$k, 1, B$k, ... , B$k, 1]  (a total of d 1's), where k ranges
# from 0 to infinity.
# To make this finite, note that we actually only need to go until d +
# k*(d-1) [which is the length of this forbidden string] is greater than n.
# So we just need k > (n-d) / (d-1).

SFP:=proc(d,n, blank:=B) local k, i:
  return {[1$d], seq([op([1, seq(op([blank$k, 1]), i=1..d-1)])], k=0..((n-d)/(d-1)+1))}:
end proc:

NdnC:=proc(d,n) local B:
  return Wn({0,1}, SFP(d,n,B), B, n):
end proc:



#############
# Problem 2 #
#############

WnX:=proc(A, FP, B, n, x) option remember:
  return WGnX(A,FP,B,n,{}, x):
end proc:


WGnX:=proc(A,FP,B,n,BFP, x) option remember: local co,a,BFP1:
  if member([],BFP union  FP) then
    RETURN(0):
  fi:
  if n=0 then
    RETURN(1):
  fi:
  co:=0:
  for a in A do
    BFP1:=AdjustBFP(a,B,FP,BFP):
    co:=co+x[a]*WGnX(A,FP,B,n-1,BFP1, x):
  od:
  return co:
end proc:


#############
# Problem 3 #
#############

adnC:=proc(d,n) local x, B:
  return degree(WnX({0,1}, SFP(d,n,B), B, n, x), x[1]):
end proc:




#############
# Problem 4 #
#############

# The sequence [seq(Wn(A,FP,B,n), n=1..infinity)] is C-finite because it is
# the linear combination of sequences of the form  WGn(A, FP, B, n, BFP).
# Moreover, we claim each such WGn is C-finite, which will complete the proof.
#
# To see that each WGn is C-finite, note that by conditioning on at most
# the first   max |b|, for b in BFP   elements of a word, we can write WGn(...) as a linear combination of the sum
# of at most this many terms of the form WGn(A, FP, B, n2, BFP2), where n2<n,
# and this combination is independent of n (provided n is larger than some constant)
#
# Moreover, the number of such terms needed in the sum to express WGn(A, ... )
# is finite and bounded by   max |f| for f in FP.
# Furthermore, the number of subsets, BFP2, that are possible to appear in these forms is also finite and bounded in terms of FP and A, which is constant.
# 
# Thus, by finiteness, these equations for WGn(A, FP, B, n, BFP) as BFP varies
# can be used to express
# Wn(A, FP, B, n) as n varies, and when n is larger than some constant
# depending on FP and A, the terms (Wn(A, FP, B, n), n=n0..infinity) will be
# annihilated by some C-finite sort of difference operator, since each of the
# terms WnG(A, FP, B, n, BFP) will be.



#############
# Problem 5 #
#############

WnR:=proc(A, FP, B, N) local n:
  return GuessRec([seq(Wn(A, FP, B, n), n=1..N)]):
end proc:


# WnR({1,2},{[1,1]},B,60);
# outputs   [[-1,-1], [2,3]]

# WnR({0,1},{[1,1,1],[0,0,0]},B,60);
# outputs [[-1,-1], [2,4]]

# WnR({0,1},{[1,1,1,1],[0,0,0,0]},B,60);
# outputs [[-1,-1,-1], [2,4,8]]

# WnR({a,b,c,d,e,f,g,h,i,j,k,l,m,n,o,p,q,r,s,t,u,w,x,y,z}, {[d,o,r,o,n],[d,B,o,B,r,B,o,B,n]},B,100);
# outputs [[-25, 0, 0, 0, 2, -25, 625, -15625, 390625, -1, 25, 0, 0, 0, -3, 25, -1250, 0, 0, 0, 0, 0, 0, 0, 1],
#          [25, 625, 15625, 390625, 9765624, 244140575, 6103513750, 152587828125, 3814694921875, 95367353515627, 2384183349609500,
#           59604571533209375, 1490113983154546875, 37252841949473046875, 931320858002089843747,
#           23283016681684814452875, 582075297832952880844375, 14551879465595245360609375, 363796912134181976284375000,
#           9094920940712451933369140630, 227372976951768398226562500500, 5684323259643375871356201207500,
#           142108052387319505134277345703125, 3552700582089014348545074564453125, 88817496362379752003192906376953117]]




################
# From Class 6 #
################

#GuessRec1(L,r): inputs a sequence of numbers, L, and a pos.
#integer r, and tries to guess a linear recurrence equation
#of order r satisfied by it
#L(n)+c1*L(n-1)+...+cr*L(n-r)=0 for all n within the range
#The  output would be given  as a list of r numbers
#[[c1,c2,...,cr],InitialConditions]
#For example, the Fibonacci sequence is
#[[-1,-1], [1,1]] 
GuessRec1:=proc(L,r) local c,i,var,eq:

if nops(L)<=2*r+6 then
 RETURN(FAIL):
fi:


var:={seq(c[i],i=1..r)}:


eq:={seq( L[n]+add(c[i]*L[n-i],i=1..r)=0 , n=r+1..nops(L))}:

var:=solve(eq,var):

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

[[seq(subs(var, c[i]),i=1..r)], L[1..r]   ] ;

end:


#GuessRec(L): inputs a sequence of numbers, L, 
# and tries to guess a linear recurrence equation
#of some order, r, satisfied by it
#L(n)+c1*L(n-1)+...+cr*L(n-r)=0 for all n within the range
#The  output would be given  as a list of r numbers
#[[c1,c2,...,cr],InitialConditions]
#For example, the Fibonacci sequence is
#[[-1,-1], [1,1]] 
GuessRec:=proc(L) local r,ans:

for r from 1 to (nops(L)-6)/2 do
 ans:=GuessRec1(L,r):
 if ans<>FAIL then
  RETURN(ans):
 fi:
od:

FAIL:
end:



#################
# From Class 19 #
#################

#Dedicated to Endre Szemeredi on the occaison of the Abel Prize 2012

with(combinat):
#Let a_d(n) be the largest size of a subset of {1, ...,n}
#avoiding arthmetical pogressions of length d
#(Szemeredi's famous theorem says that a_d(n)/n->0)
#d=3 is it the already deep Roth theorem who proved that
#a_3(n)/n<=C/log(log(n))




#Adn(d,n): Inputs positive integers d and n and outputs
#the set of all subsets of {1, ...,n} that DO NOT have
#an arithmetical progression of length d
Adn:=proc(d,n) local S,s1,i,T:

S:=powerset(n):
T:={}:

for s1 in S do 

if not IsBad(d,s1) then

#T:=T union {s1}:
 
T:={op(T), s1}:

fi:
 
od:

T:


end:




#IsBad(d,S): Does the set of pos. integers  S contain
#an AP of length d?
IsBad:=proc(d,S) local m,S1,dif,PTM,j:
option remember:

if nops(S)<d then
RETURN(false):
fi:

m:=max(S):

S1:=S minus {m}:

if IsBad(d,S1) then
  RETURN(true):
fi:

for dif from 1 to trunc(m/(d-1)) do
PTM:={seq(m-j*dif, j=1..d-1)}:

if PTM minus S1={} then
  RETURN(true):
fi:

od:


false:
end:

#adnS(d,n): a stupid way of computing a_d(n) above
adnS:=proc(d,n) local Y,y:
Y:=Adn(d,n):

max(seq(nops(y),y in Y)):

end:


#NdnS(d,n): a stupid way of computing the number of 
#subsets of {1, ...,n} that do not contain an AP
#of length d
NdnS:=proc(d,n) local Y,y:
Y:=Adn(d,n):

nops(Y):

end:




#Wn(A,FP,B,n): inputs
#(i) A: a set of symbols (the alphabet)
#(ii) FP: a set of forbidden patterns, given as a list
#where B denotes Blank. For example, if FP contains
#[1,B,1,B,1] then our words may not contain
#10101, 10111, 11101, 11111
#(iii) A letter B
#n a pos. integer
#Outputs:
#The number of words in the alphabet A of length n that
#AVOID all the patterns of FP.
#For example, 
#Wn({0,1},{[1,B,1]},B,3); should return
#6
Wn:=proc(A,FP,B,n) 
option remember:
WGn(A,FP,B,n,{});
end:


#




#AdnM(d,n): mediocre version of Adn(d,n)
#Inputs positive integers d and n and outputs
#the set of all subsets of {1, ...,n} that DO NOT have
#an arithmetical progression of length d
AdnM:=proc(d,n) local S,T,t1,s1:
option remember:

if n=1 then 
RETURN({{},{1}}):
fi:

S:=AdnM(d,n-1):

T:=S:

for t1 in T do 

s1:=t1 union {n}:

if not IsBad(d,s1) then

#S:=S union {s1}:
 
S:={op(S), s1}:

fi:
 
od:

S:


end:

#adnM(d,n): a mediocre way of computing a_d(n) above
adnM:=proc(d,n) local Y,y:
Y:=AdnM(d,n):

max(seq(nops(y),y in Y)):

end:


#NdnM(d,n): a stupid way of computing the number of 
#subsets of {1, ...,n} that do not contain an AP
#of length d
NdnM:=proc(d,n) local Y,y:
Y:=AdnM(d,n):

nops(Y):

end:

#WGn(A,FP,B,n,BFP): inputs
#(i) A: a set of symbols (the alphabet)
#(ii) FP: a set of forbidden patterns, given as a list
#where B denotes Blank. For example, if FP contains
#[1,B,1,B,1] then our words may not contain
#10101, 10111, 11101, 11111
#(iii) A letter B
#(iv) n a pos. integer
#(v): Another set of patterns (same format at for FP)
#Outputs:
#The number of words in the alphabet A of length n that
#AVOID all the patterns of FP and IN ADDITION AVOID
#the set of patterns BFP at the beginning
#For example, 
#Wn({0,1},{[1,B,1]},B,3); should return
#6
WGn:=proc(A,FP,B,n,BFP) local co,a,BFP1:
option remember:
if member([],BFP union  FP) then
  RETURN(0):
fi:

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


co:=0:

for a in A do

BFP1:=AdjustBFP(a,B,FP,BFP):

co:=co+WGn(A,FP,B,n-1,BFP1):

od:

co:

end:


#AdjustBFP(a,B,FP,BFP): inputs a letter a, a symbol B for
#blank, a set of forbidden patterns FP, and another set
#of forbidden patterns at the beginning
#outputs the new BFP for the chopped word assuming
#it starts with a
AdjustBFP:=proc(a,B,FP,BFP)  local  NBFP,fp:

NBFP:={}:


if member([],BFP) then
 RETURN(FAIL):
fi:

for fp in FP union BFP do

 if fp[1]=a or fp[1]=B then
   NBFP:=NBFP union {fp[2..nops(fp)]}:
  fi:

od:



NBFP:

end: