#####################################################################
## YoungT.txt Save this file as  YoungT.txt  to use it,              #
# stay in the                                                        #
## same directory, get into Maple (by typing: maple <Enter> )        #
## and then type:  read YountT.txt <Enter>                           #
## Then follow the instructions given there                          #
##                                                                   #
## Written by Doron Zeilberger, Rutgers University ,                 #
##  DoronZeil at gmail dot com                                       # 
######################################################################

with(plots):
with(combinat):
Digits:=40:

print(`First Written May  2020: tested for Maple 18 `):
print(`Version  of May 2020  `):
print():
print(`This is YoungT.txt, one of the Maple packages that `):
print(`accompany Maneul Kauers's and Doron Zeilberger's  article: `):
print(` Counting Standard Young Tableaux with Restricted Runs`):
print(`Available from:`):
print(` http://sites.math.rutgers.edu/~zeilberg/mamarim/mamarimhtml/cyt.html .`):
print():
print(`The most current version is available on WWW at:`):
print(` http://sites.math.rutgers.edu/~zeilberg/tokhniot/YoungT.txt .`):
print(`Please report all bugs to: DoronZeil at gmail dot com .`):



print():
print(`For a list of the  STORY functions,`):
print(` type "ezraStory();". For specific help type "ezra(procedure_name);" `):
print():


print():
print(`For a list of the  CHECKING functions,`):
print(` type "ezraCheck();". For specific help type "ezra(procedure_name);" `):
print():

print():
print(`For a list of the SUPPORTING functions,`):
print(` type "ezra1();". For specific help type "ezra(procedure_name);" `):
print():


print():
print(`For a list of the procedures from FindRec.txt`):
print(` type "ezraFindRec();". For specific help type "ezraFindRec(procedure_name);" `):
print():
 

print():
print(`For general help, and a list of the MAIN functions,`):
print(` type "ezra();". For specific help type "ezra(procedure_name);" `):
#print(`For a list of the supporting functions type: ezra1();`):
print():
 




ezraCheck:=proc()
if args=NULL then

print(` The Checking procedures are :  GmResSet, MacM, Runs, YT  `):

else
ezra(args):
fi:

end:


ezraStory:=proc()
if args=NULL then

print(` The Story procedures are : GseqVA , GseqRVA `):


else
ezra(args):

fi:


end:


ezra1:=proc()
if args=NULL then

print(` The supporting procedures are : Fm, Fmi, FmiRes,FmiResR, FmRes,FmResR, FmResRViaGF, FmResViaGF, Fkx, FkxRes, FkxResR, Gm, Gmi,`):
print(`  GmiRes, GmoResR,  GmResR , MyAsy, MyAsyC, Zinn, ZinnPlus`):

else
ezra(args):

fi:


end:


ezra:=proc()
if args=NULL then

print(` YoungT.txt:  A Maple package for automatic enumeration of families of Young Tableaux (and paths) obeying certain restrictions  `):

print(`The MAIN procedures are`):
print(` Fseq, FseqR, Gseq , GseqR, InfoF, InfoFr, InfoG, InfoGr `):


elif nargs=1 and args[1]=Fseq then
print(`Fseq(RES,N): inputs a list R (let d:=nops(R)) of sets of positive integers outputs the sequence [seq(FmRes([n$d],RES), n=1..N)]. Try`):
print(`Fseq([{1},{1}],20);`);


elif nargs=1 and args[1]=Fm then
print(`Fm(m): inputs a list m of non-negative integers m, indicating a point in (let d:=nops(m)) the d-dimensional lattice`):
print(`and 1<=i<=d, outputs the number of lattice walks that end at m. Try:`):
print(`Fm([4,5,3]);`):

elif nargs=1 and args[1]=Fmi then
print(`Fmi(m,i): inputs a list m of non-negative integers m, indicating a point in (let d:=nops(m)) the d-dimensional lattice`):
print(`and 1<=i<=d, outputs the number of lattice walks that end in a run of steps parallel to the i-th axis. Try:`):
print(`Fmi([4,5,3],2);`):


elif nargs=1 and args[1]=FmiRes then
print(`FmiRes(m,i,R): inputs a list m of non-negative integers m, indicating a point in (let d:=nops(m)) the d-dimensional lattice`):
print(`and 1<=i<=d, as well as a list R of length d, of sets of positive integers`):
print(`outputs the number of lattice walks that end in a run of steps parallel to the i-th axis. and such that`):
print(`none of the runs in the i-th dimension is of length that belongs to R[i]. Try:`):
print(`FmiRes([4,5,3],2,[{1},{1},{1,2}]);`):


elif nargs=1 and args[1]=FmRes then
print(`FmRes(m,R): inputs a list m of non-negative integers m, indicating a point in (let d:=nops(m)) the d-dimensional lattice`):
print(`and 1<=i<=d, as well as a list R of length d, of sets of positive integers`):
print(`outputs the number of lattice walks that end in m`):
print(`none of the runs in the i-th dimension is of length that belongs to R[i]. Try:`):
print(`FmRes([4,5,3],[{1},{1},{1,2}]);`):


elif nargs=1 and args[1]=FmResR then
print(`FmResR(m,R,r): inputs a list m of non-negative integers m, indicating a point in (let d:=nops(m)) the d-dimensional lattice`):
print(`and 1<=i<=d, as well as a list R of length d, of affine linear expressions in the variable r, with integer coefficients.`):
print(`Outputs the number of lattice walks that end in m`):
print(`none of the runs in the i-th dimension is of length that belongs to the range of R[i]. Try:`):
print(`FmResR([4,5,3],[2*r+2,2*r+1,2*r+2],r);`):



elif nargs=1 and args[1]=FmResRViaGF then
print(`FmResRViaGF(m,R,r): like FmResR(m,R,r) (q.v.) but via the generating function FkxResR(k,R,r,x) (q.v.). Try:`):
print(`FmResRViaGF([3,5], [2*r+1,3*r+1],r);`):

elif nargs=1 and args[1]=FmResViaGF then
print(`FmResViaGF(m,R): like FmRes(m,R) (q.v.) but via the generating function FkxRes(k,R,x) (q.v.). Try:`):
print(`FmResViaGF([4,5,3],[{1},{1},{1,2}]);`):

elif nargs=1 and args[1]=Fkx then
print(`Fkx(k,x): The generating function of all words in {1,...,k}. Try:`):
print(`Fkx(k,x);`):

elif nargs=1 and args[1]=FkxRes then
print(`FkxRes(k,R,x): The generating function of all words in {1,...,k} with no i-runs of length in the set R[i].`):
print(`Try:`):
print(`FkxRes(2,[{1},{2}],x);`):


elif nargs=1 and args[1]=FkxResR then
print(`FkxResR(k,R,r,x): The generating function of all words in {1,...,k} with no i-runs of length in the range of R[i].`):
print(`where R[i] is an affine linear expression in r.`):
print(`Try:`):
print(`FkxResR(2,[2*r+1,3*r+1],r,x);`):

elif nargs=1 and args[1]=FseqR then
print(`FseqR(R,r,N): inputs a list R (let d:=nops(R)) of affine linear expressions in the variable r,`):
print(` outputs the sequence [seq(FmResR([n$d],R,r), n=1..N)]. Try:`):
print(`FseqR([2*r+2,2*r+2,2*r+2],r,30);`):

elif nargs=1 and args[1]=Gm then
print(`Gm(m): inputs a list m of non-negative integers m, indicating a point in (let d:=nops(m)) the d-dimensional lattice`):
print(`in the region x[1]>=...>=x[d]>=0 `):
print(`and 1<=i<=d, outputs the number of lattice walks that end at m. Try:`):
print(`Gm([5,3,3]);`):



elif nargs=1 and args[1]=Gmi then
print(`Gmi(m,i): inputs a list m of non-negative integers m, indicating a point in (let d:=nops(m)) the d-dimensional lattice`):
print(`in the region x[1]>=...>=x[d]>=0 `):
print(`and 1<=i<=d, outputs the number of lattice walks that end in a run of steps parallel to the i-th axis. Try:`):
print(`Gmi([5,4,3],2);`):


elif nargs=1 and args[1]=GmiRes then
print(`GmiRes(m,i,R): inputs a list m of non-negative integers m, indicating a point in (let d:=nops(m)) the d-dimensional lattice`):
print(`and 1<=i<=d, as well as a list R of length d, of sets of positive integers`):
print(`outputs the number of lattice walks that end in a run of steps parallel to the i-th axis. `):
print(`and stay in x[1]>=..>=x[d]>=0`):
prin(`and such that`):
print(`none of the runs in the i-th dimension is of length that belongs to R[i]. Try:`):
print(`FmiRes([4,5,3],2,[{1},{1},{1,2}]);`):


elif nargs=1 and args[1]=GmiResR then
print(`GmiResR(m,i,R,r): inputs a list m of non-negative integers m, indicating a point in (let d:=nops(m)) the d-dimensional lattice`):
print(`and 1<=i<=d, as well as a list R of length d, of affine-linear expressions in r`):
print(`outputs the number of lattice walks that end in a run of steps parallel to the i-th axis. `):
print(`and stay in x[1]>=..>=x[d]>=0`):
prin(`and such that`):
print(`none of the runs in the i-th dimension is of length that belongs to the range of R[i]. Try:`):
print(`FmiResR([4,5,3],2,[2*r+1,3*r+1,2*r],r);`):


elif nargs=1 and args[1]=GmRes then
print(`GmRes(m,R): inputs a list m of non-negative integers m, indicating a point in (let d:=nops(m)) the d-dimensional lattice`):
print(`and 1<=i<=d, as well as a list R of length d, of sets of positive integers`):
print(`outputs the number of lattice walks that end in m`):
print(`and stay in x[1]>=..>=x[d]>=0`):
print(`none of the runs in the i-th dimension is of length that belongs to R[i]. Try:`):
print(`GmRes([5,3,3],[{1},{1},{1,2}]);`):


elif nargs=1 and args[1]=GmResR then
print(`GmResR(m,R,r): inputs a list m of non-negative integers m, indicating a point in (let d:=nops(m)) the d-dimensional lattice`):
print(`and 1<=i<=d, as well as a list R of length d, of affine linear expressions in r`):
print(`outputs the number of lattice walks that end in m`):
print(`and stay in x[1]>=..>=x[d]>=0`):
print(`none of the runs in the i-th dimension is of length that belongs to the range of R[i]. Try:`):
print(`GmResR([5,3,3],[2*r+1,2*r+1,2*r+1],r);`):



elif nargs=1 and args[1]=GmResSet then
print(`GmResSet(m,R): The set of Young Tableau of shaper m such that the i-th row has no run-lenghths in R[i]. Try:`):
print(`GmResSet([4,3,2],[{1},{1},{1}]);`):

elif nargs=1 and args[1]=Gseq then
print(`Gseq(R,N): inputs a list R (let d:=nops(R)) of sets of positive integers outputs the sequence [seq(GmRes([n$d],R), n=1..N)]. Try`):
print(`Gseq([{1},{1}],20);`);

elif nargs=1 and args[1]=GseqVA then
print(`GseqVA(R,N,n): Verbose form of Gseq(R,N) (q.v.) and also gives the asympototics to order 2. N has to be at least 50`):
print(`Try: `):
print(`GseqVA([{1},{1},{1}],50,n);`):

elif nargs=1 and args[1]=GseqR then
print(`GseqR(R,r,N): inputs a list R (let d:=nops(R)) of sets of positive integers outputs the sequence [seq(FmRes([n$d],R), n=1..N)]. Try`):
print(`GseqR([2*r+1,2*r+1],r,20);`);


elif nargs=1 and args[1]=GseqRVA then
print(`GseqRVA(R,r,N,n): Verbose form of GseqR(R,r,N) (q.v.) and also gives the asympototics to order 2. N has to be at least 50`):
print(`Try: `):
print(`GseqRVA([2*r+2,2*r+2,2*r+2],r,50,n);`):




elif nargs=1 and args[1]=InfoF then
print(`InfoF(RES,N0,n,a,MaxC,N1): Inputs`):
print(`(i) a list RES of sets of posivive integers, let k=nops(RES). outputs`):
print(`(ii) a fairly large pos. integer N0`):
print(`(iii) symbols n and a`):
print(`(iv) a pos. integer MaxC`):
print(`(v): A large positive integer N1`):
print(`Outputs a list consisting of`):
print(`(i) The first N0 terms of the sequence "number of walks in k-dim hyper-cubic lattice`):
print(` from the origin to [n,...,n] without runs in the i-th direction of length in R[i] `):
print(` (ii) A linear recurrence equation with poly coefficients of complexity <=MaxC `):
print(`satisfied by the sequence in terms of a(n), followed by the initial conditions . If none found it only returns the above list`):
print(`(iii) The value of the N1-th term of the sequence. Try`):
print(`InfoF([{1},{1}],50,n,a,10,100);`):


elif nargs=1 and args[1]=InfoFr then
print(`InfoFr(RES,r,N0,n,a,MaxC,N1): Inputs`):
print(`(i) a list RES of affine linear expressions in r, let k=nops(RES). outputs`):
print(`(ii) a fairly large pos. integer N0`):
print(`(iii) symbols n and a`):
print(`(iv) a pos. integer MaxC`):
print(`(v): A large positive integer N1`):
print(`Outputs a list consisting of`):
print(`(i) The first N0 terms of the sequence "number of walks in k-dim hyper-cubic lattice`):
print(` from the origin to [n,...,n] without runs in the i-th direction of length in the range of R[i] `):
print(` (ii) A linear recurrence equation with poly coefficients of complexity <=MaxC `):
print(`satisfied by the sequence in terms of a(n), followed by the initial conditions . If none found it only returns the above list`):
print(`(iii) The value of the N1-th term of the sequence. Try`):
print(`InfoFr([r+2,r+2],r,50,n,a,10,100);`):



elif nargs=1 and args[1]=InfoG then
print(`InfoG(RES,N0,n,a,MaxC,N1): Inputs`):
print(`(i) a list RES of sets of posivive integers, let k=nops(RES). outputs`):
print(`(ii) a fairly large pos. integer N0`):
print(`(iii) symbols n and a`):
print(`(iv) a pos. integer MaxC`):
print(`(v): A large positive integer N1`):
print(`Outputs a list consisting of`):
print(`(i) The first N0 terms of the sequence "number of walks in k-dim hyper-cubic lattice`):
print(` from the origin to [n,...,n] staying in x[1]>=...>=x[k] without runs in the i-th direction of length in R[i] `):
print(` (ii) A linear recurrence equation with poly coefficients of complexity <=MaxC `):
print(`satisfied by the sequence in terms of a(n), followed by the initial conditions . If none found it only returns the above list`):
print(`(iii) The value of the N1-th term of the sequence. Try`):
print(`InfoG([{1},{1}],50,n,a,10,100);`):


elif nargs=1 and args[1]=InfoGr then
print(`InfoGr(RES,r,N0,n,a,MaxC,N1): Inputs`):
print(`(i) a list RES of affine linear expressions in r, let k=nops(RES). outputs`):
print(`(ii) a fairly large pos. integer N0`):
print(`(iii) symbols n and a`):
print(`(iv) a pos. integer MaxC`):
print(`(v): A large positive integer N1`):
print(`Outputs a list consisting of`):
print(`(i) The first N0 terms of the sequence "number of walks in k-dim hyper-cubic lattice`):
print(` from the origin to [n,...,n] staying in x[1]>=...>=x[k], without runs in the i-th direction of length in the range of R[i] `):
print(` (ii) A linear recurrence equation with poly coefficients of complexity <=MaxC `):
print(`satisfied by the sequence in terms of a(n), followed by the initial conditions . If none found it only returns the above list`):
print(`(iii) The value of the N1-th term of the sequence. Try`):
print(`InfoGr([r+2,r+2],r,50,n,a,10,100);`):

elif nargs=1 and args[1]=MacM then
print(`MacM(a): MacMahon's formula for the number of Young tableaux of shape a. Try:`):
print(`MacM([3,3,2]);`):

elif nargs=1 and args[1]=MyAsy then
print(`MyAsy(L,n,k): tries to fit the sequence L to be of the form  C*mu^n*n^theta*(1+c1/n+...+ck/n^k). Try:`):
print(`L:=Gseq([{1},{1}],50); MyAsy(L,n,3);`):


elif nargs=1 and args[1]=MyAsyC then
print(`MyAsyC(L,n,k,theta,mu): tries to fit the sequence L to be of the form  C*mu^n*n^theta*(1+c1/n+...+ck/n^k). `):
print(`when mu and theta are given. Try:`):
print(`L:=Gseq([{},{}],50); MyAsyC(L,n,3,-3/2,4);`):

elif nargs=1 and args[1]=Runs1 then
print(`Runs1(L): Given  an incerasing list of postive integers outputs the sequence of runs. Try:`):
print(`Runs1([5,7,8,10]);`):

elif nargs=1 and args[1]=YT then
print(`YT(L): The set of Young Tableaux of shape L. Try:`):
print(`YT([2,2,2]);`):


elif nargs=1 and args[1]=Zinn then
print(`Zinn(L): estimates [theta,mu] such that a(n) is asymptotic to mu^n*n^theta`):
print(`Try:`):
print(`Zinn(SeqG([{1},{1},{1}],30);`):

elif nargs=1 and args[1]=ZinnPlus then
print(`ZinnPlus(resh,n): An estimate for the asympotics of resh[n] in the form C*mu^n*n^theta. Try:`):
print(`L:=Fseq([{},{}],100): ZinnPlus(L,n);`):

else
 print(`There is no such thing as`, args):

fi:


end:


   
 
##start from Findrec
ezraFindRec:=proc()
if args=NULL then

print(` FindRec.txt: A Maple package for empirically guessing partial recurrence`):
print(`equations satisfied by Discrete Functions of ONE Variable`):
print():
print(`For help with a specific procedure, type "ezra(procedure_name);"`):
print(`Contains procedures:  `):
print(`  findrec, Findrec, FindrecF, Zinn, ZinnPlus `):
print():

elif nargs=1 and args[1]=findrec then
print(`findrec(f,DEGREE,ORDER,n,N): guesses a recurrence operator annihilating`):
print(`the sequence f of degree DEGREE and order ORDER.`):
print(`For example, try: findrec([seq(i,i=1..10)],0,2,n,N);`):

elif nargs=1 and args[1]=Findrec then
print(`Findrec(f,n,N,MaxC): Given a list f tries to find a linear recurrence equation with`):
print(`poly coffs. ope(n,N), where n is the discrete variable, and N is the shift operator `):
print(`of maximum DEGREE+ORDER<=MaxC`):
print(`e.g. try Findrec([1,1,2,3,5,8,13,21,34,55,89],n,N,2);`):

elif nargs=1 and args[1]=FindrecF then
print(`FindrecF(f,n,N): Given a function f of a single variable tries to find a linear recurrence equation with`):
print(`poly coffs. .g. try FindrecF(i->i!,n,N);`):


elif nargs=1 and args[1]=SeqFromRec then
print(`SeqFromRec(ope,n,N,Ini,K): Given the first L-1`):
print(`terms of the sequence Ini=[f(1), ..., f(L-1)]`):
print(`satisfied by the recurrence ope(n,N)f(n)=0`):
print(`extends it to the first K values`):
print(`For example, try:`):
print(`SeqFromRec(N-n-1,n,N,[1],10);`):


fi:
 
end:


Zinn:=proc(resh)
local s1,s2,n:
n:=nops(resh)-2:
s1:=sn(resh,n):
s2:=sn(resh,n-1):
evalf(2*(s1+s2)/(s1-s2)^2),
evalf(sqrt(op(n+1,resh)/op(n-1,resh))*exp(-(s1+s2)/((s1-s2)*s1))):
end:
 
 
sn:=proc(resh,n):
-1/log(op(n+1,resh)*op(n-1,resh)/op(n,resh)^2):
#evalf(%):
end:
 

###Findrec 
#findrec(f,DEGREE,ORDER,n,N): guesses a recurrence operator annihilating
#the sequence f of degree DEGREE and order ORDER
#For example, try: findrec([seq(i,i=1..10)],0,2,n,N);
findrecVerbose:=proc(f,DEGREE,ORDER,n,N)
local ope,var,eq,i,j,n0,kv,var1,eq1,mu,a:
if (1+DEGREE)*(1+ORDER)+3+ORDER>nops(f) then
#ERROR(`Insufficient data for a recurrence of order`,ORDER, `degree`,DEGREE):
RETURN(FAIL):
fi:
ope:=0:
var:={}:
 
for i from 0 to ORDER do
 for j from 0 to DEGREE do
  ope:=ope+a[i,j]*n^j*N^i:
  var:=var union {a[i,j]}:
 od:
od:
    
eq:={}:
 
for n0 from 1 to (1+DEGREE)*(1+ORDER)+2 do
  eq1:=0:
 
  for i from 0 to ORDER do
     eq1:=eq1+subs(n=n0,coeff(ope,N,i))*op(n0+i,f):
  od:
 
   eq:= eq union {eq1}:
od:
 
var1:=solve(eq,var):
 
kv:={}:
 
for i from 1 to nops(var1) do
 mu:=op(i,var1):
 
if op(1,mu)=op(2,mu) then
   kv:= kv union {op(1,mu)}:
 fi:
 
od:

ope:=subs(var1,ope):

if ope=0 then
  RETURN(FAIL):
fi:

ope:={seq(coeff(expand(ope),kv[i],1),i=1..nops(kv))} minus {0}:

if nops(ope)>1 then
print(`There is some slack, there are `, nops(ope)):
print(ope):
RETURN(Yafe(ope[1],N)[2]):
elif nops(ope)=1 then
RETURN(Yafe(ope[1],N)[2]):
else
 RETURN(FAIL):
fi:

end:
 
#findrec(f,DEGREE,ORDER,n,N): guesses a recurrence operator annihilating
#the sequence f of degree DEGREE and order ORDER
#For example, try: findrec([seq(i,i=1..10)],0,2,n,N);
findrec:=proc(f,DEGREE,ORDER,n,N)
local ope,var,eq,i,j,n0,kv,var1,eq1,mu,a,i1,n1:
option remember:


if not findrecEx(f,DEGREE,ORDER,ithprime(20)) then
 RETURN(FAIL):
fi:

if not findrecEx(f,DEGREE,ORDER,ithprime(40)) then
 RETURN(FAIL):
fi:

if not findrecEx(f,DEGREE,ORDER,ithprime(80)) then
 RETURN(FAIL):
fi:


if (1+DEGREE)*(1+ORDER)+5+ORDER>nops(f) then
#ERROR(`Insufficient data for a recurrence of order`,ORDER, `degree`,DEGREE):
RETURN(FAIL):
fi:
ope:=0:
var:={}:
 
for i from 0 to ORDER do
 for j from 0 to DEGREE do
  ope:=ope+a[i,j]*n^j*N^i:
  var:=var union {a[i,j]}:
 od:
od:
    
eq:={}:
 
for n0 from 1 to (1+DEGREE)*(1+ORDER)+4 do
  eq1:=0:
 
  for i from 0 to ORDER do
     eq1:=eq1+subs(n=n0,coeff(ope,N,i))*op(n0+i,f):
  od:
 
   eq:= eq union {eq1}:
od:
 
var1:=solve(eq,var):
 
kv:={}:
 
for i from 1 to nops(var1) do
 mu:=op(i,var1):
 
if op(1,mu)=op(2,mu) then
   kv:= kv union {op(1,mu)}:
 fi:
 
od:

ope:=subs(var1,ope):

if ope=0 then
  RETURN(FAIL):
fi:

ope:={seq(coeff(expand(ope),kv[i],1),i=1..nops(kv))} minus {0}:

if nops(ope)>1 then
RETURN(Yafe(ope[1],N)[2]):
elif nops(ope)=1 then
ope:=ope[1]:

if {seq(add(subs(n=n1+i1,coeff(ope,N,i1))*f[n1+i1],i1=0..degree(ope,N)),n1=1..nops(f)-degree(ope,N))}<>{0} then
 RETURN(FAIL):
fi:

 ope:=Yafe(ope,N)[2]:


if SeqFromRec(ope,n,N,[op(1..degree(ope,N),f)],nops(f))<>f then
 RETURN(FAIL):
else
RETURN(ope):
fi:
else
 RETURN(FAIL):
fi:

end:
 


Yafe:=proc(ope,N) local i,ope1,coe1,L: 
if ope=0 then
 RETURN(1,0):
fi:
ope1:=expand(ope):
L:=degree(ope1,N):
coe1:=coeff(ope1,N,L):
ope1:=normal(ope1/coe1):
ope1:=normal(ope1):
ope1:=
convert(
[seq(factor(coeff(ope1,N,i))*N^i,i=ldegree(ope1,N)..degree(ope1,N))],`+`):
factor(coe1),ope1:
end:


#Findrec(f,n,N,MaxC): Given a list f tries to find a linear recurrence equation with
#poly coffs.
#of maximum DEGREE+ORDER<=MaxC
#e.g. try Findrec([1,1,2,3,5,8,13,21,34,55,89],n,N,2);
Findrec:=proc(f,n,N,MaxC)
local DEGREE, ORDER,ope,L:

for L from 0 to MaxC do
for ORDER from 0 to L do
 DEGREE:=L-ORDER:
if (2+DEGREE)*(1+ORDER)+4>=nops(f) then
# print(`Insufficient data for degree`, DEGREE, `and order `,ORDER):
 RETURN(FAIL):
fi:
 ope:=findrec([op(1..(2+DEGREE)*(1+ORDER)+4,f)],DEGREE,ORDER,n,N):
     if ope<>FAIL and degree(ope,N)<>0 then
       RETURN(ope):
      fi:
 od:
od:
FAIL:

end:





 
#SeqFromRec(ope,n,N,Ini,K): Given the first L-1
#terms of the sequence Ini=[f(1), ..., f(L-1)]
#satisfied by the recurrence ope(n,N)f(n)=0
#extends it to the first K values
SeqFromRec:=proc(ope,n,N,Ini,K)
local ope1,gu,L,n1,j1,ka,i1:
ope1:=ope:

L:=degree(ope1,N):
if nops(Ini)<>L then
 ERROR(`Ini should be of length`, L):
fi:

ka:={seq(subs(n=i1,lcoeff(ope1,N)),i1=1..K)}:

if member(0,ka) then
 RETURN(FAIL):
fi:

ope1:=expand(subs(n=n-L,ope1)/N^L):

gu:=Ini:

for n1 from nops(Ini)+1 to K do
if member(0, { seq( subs(n=n1,denom(coeff(ope1,N,-j1)  ))  ,j1=1..L )}) then
 RETURN(FAIL):
else
gu:=[op(gu), -add(gu[nops(gu)+1-j1]*subs(n=n1,coeff(ope1,N,-j1)),
j1=1..L)]:
fi:
od:

gu:

end:
 
#end Findrec


with(linalg):

#findrecEx(f,DEGREE,ORDER,m1): Explores whether thre
#is a good chance that there is a recurrence of degree DEGREE
#and order ORDER, using the prime m1
#For example, try: findrecEx([seq(i,i=1..10)],0,2,n,N,1003);
findrecEx:=proc(f,DEGREE,ORDER,m1)
local ope,var,eq,i,j,n0,eq1,a,A1,
D1,E1,Eq,Var,f1,n,N:
option remember:
f1:=f mod m1:
if (1+DEGREE)*(1+ORDER)+5+ORDER>nops(f) then
#ERROR(`Insufficient data for a recurrence of order`,ORDER, `degree`,DEGREE):
RETURN(FAIL):
fi:
ope:=0:
var:={}:
 
for i from 0 to ORDER do
 for j from 0 to DEGREE do
  ope:=ope+a[i,j]*n^j*N^i:
  var:=var union {a[i,j]}:
 od:
od:
    
eq:={}:
 
for n0 from 1 to (1+DEGREE)*(1+ORDER)+4 do
  eq1:=0:
 
  for i from 0 to ORDER do
     eq1:=eq1+subs(n=n0,coeff(ope,N,i))*op(n0+i,f1) mod m1:
  od:
 
   eq:= eq union {eq1}:
od:


Eq:= convert(eq,list):
Var:= convert(var,list):


D1:=nops(Var):
E1:=nops(Eq):
if E1<D1 then
  RETURN(true):
fi:

A1:=matrix(D1,D1):

for i from 1 to D1-1 do
for j from 1 to D1 do
  A1[i,j]:=coeff(Eq[i],Var[j]):
od:
od:

for j from 1 to nops(Var) do
  A1[D1,j]:=coeff(Eq[D1],Var[j]):
od:

if det(A1) mod m1 <>0  then
 RETURN(false):
fi:

if E1-D1>=1 then
 for j from 1 to nops(Var) do
  A1[D1,j]:=coeff(Eq[D1+1],Var[j]):
 od:

if det(A1) mod m1 <>0  then
 RETURN(false):
fi:
fi:

if E1-D1>=2 then
 for j from 1 to nops(Var) do
  A1[D1,j]:=coeff(Eq[D1+2],Var[j]):
 od:

if det(A1) mod m1 <>0  then
 RETURN(false):
fi:
fi:

true:

end:
###End from Findrec



#OpeToRec(ope,n,N,INI,a): inputs a recurrence operators ope in n and N and initial coditions, write it nicel
# in terms of a(n)=Combination of a(n-i), and the initial conditions. Try:
#OpeToRec(N^2-N-1,n,N,[1,1],a);
OpeToRec:=proc(ope,n,N,INI,a) local ORD,i,ope1:

ORD:=degree(ope,N):

ope1:=expand(subs(n=n-ORD,ope)/N^ORD):

ope1:=1-expand(ope1/coeff(ope1,N,0)):

[a(n)=add(factor(coeff(ope1,N,-i))*a(n-i),i=1..ORD), [seq(a(i)=INI[i],i=1..ORD)]]:

end:


#Fmi(m,i): inputs a list m of non-negative integers m, indicating a point in (let d:=nops(m)) the d-dimensional lattice
#and 1<=i<=d, outputs the number of lattice walks that end in a run of steps parallel to the i-th axis.
#Try:
#Fmi([4,5,3],2);
Fmi:=proc(m,i) local d,j,k,i1:
option remember:
d:=nops(m):

if not (1<=i and i<=d) then
 RETURN(FAIL):
fi:


if {seq(m[i1],i1=1..i-1),seq(m[i1],i1=i+1..d)}={0} then
 RETURN(1):
fi:

add(add(Fmi([op(1..i-1,m), j,op(i+1..d,m)],k),j=0..m[i]-1),k=1..i-1)+
add(add(Fmi([op(1..i-1,m), j,op(i+1..d,m)],k),j=0..m[i]-1),k=i+1..d):
end:

#Fm(m): inputs a list m of non-negative integers m, indicating a point in (let d:=nops(m)) the d-dimensional lattice
# outputs the number of positive lattice walks from the origin to m. Try:
#Fm([4,5,3]);
Fm:=proc(m) local i:
option remember:
add(Fmi(m,i),i=1..nops(m)):
end:



#Gmi(m,i): inputs a list m of non-negative integers m, indicating a point in (let d:=nops(m)) the d-dimensional lattice
#and the regin m[1]>=...>=m[d]>=0
#and 1<=i<=d, outputs the number of lattice walks that end in a run of steps parallel to the i-th axis.
#Try:
#Gmi([4,5,3],2);
Gmi:=proc(m,i) local d,j,k,i1:
option remember:
d:=nops(m):

if not (1<=i and i<=d) then
 RETURN(FAIL):
fi:

if min(seq(m[i1]-m[i1+1],i1=1..d-1))<0 then
 RETURN(0):
fi:


if i=1 and {seq(m[i1],i1=2..d)}={0} then
 RETURN(1):
fi:

add(add(Gmi([op(1..i-1,m), j,op(i+1..d,m)],k),j=0..m[i]-1),k=1..i-1)+
add(add(Gmi([op(1..i-1,m), j,op(i+1..d,m)],k),j=0..m[i]-1),k=i+1..d):
end:


#Gm(m): inputs a list m of non-negative integers m, indicating a point in (let d:=nops(m)) the d-dimensional lattice
# outputs the number of positive lattice walks from the origin to m. Try:
#Gm([4,5,3]);
Gm:=proc(m) local i:
option remember:
add(Gmi(m,i),i=1..nops(m)):
end:


#FmiRes(m,i,R): inputs a list m of non-negative integers m, indicating a point in (let d:=nops(m)) the d-dimensional lattice
#and 1<=i<=d,  and also a list of length d, of sets of positive integers,
#outputs the number of lattice walks that end in a run of steps parallel to the i-th axis.
#and that each consecutive run in the i-th axis can't belong to R[i]. Try
#Try:
#FmiRes([5,2,3],2,[{1},{1,2},{3}] );
FmiRes:=proc(m,i,R) local i1,d,j,S,s,k:
option remember:
d:=nops(m):

if nops(R)<>d then
 RETURN(FAIL):
fi:

if not (1<=i and i<=d) then
 RETURN(FAIL):
fi:


if ({seq(m[i1],i1=1..i-1),seq(m[i1],i1=i+1..d)}={0} or {seq(m[i1],i1=1..i-1),seq(m[i1],i1=i+1..d)}={}) then
if member(m[i],R[i]) then
 RETURN(0):
else
 RETURN(1):
fi:

fi:

S:={seq(j,j=1..m[i])}  minus R[i]:

add(add(FmiRes([op(1..i-1,m), m[i]-s,op(i+1..d,m)],k,R),s in S),k=1..i-1)+
add(add(FmiRes([op(1..i-1,m), m[i]-s,op(i+1..d,m)],k,R), s in S),k=i+1..d):
end:


#FmRes(m,R): inputs a list m of non-negative integers m, indicating a point in (let d:=nops(m)) the d-dimensional lattice
#and a list of sets R of length d
# outputs the number of positive lattice walks from the origin to m such that no runs in the x[i] direction can
#have a length that belongs to R[i]. Try:
#FmRes([4,5,3],[{1},{2},{1}]);
FmRes:=proc(m,R) local i:
option remember:
add(FmiRes(m,i,R),i=1..nops(m)):
end:





#GmiRes(m,i,R): inputs a list m of non-negative integers m, indicating a point in (let d:=nops(m)) the d-dimensional lattice
#and 1<=i<=d,  and also a list of length d, of sets of positive integers,
#outputs the number of lattice walks  in x[1]>=..x[d]>=0
#that end in a run of steps parallel to the i-th axis.
#and that each consecutive run in the i-th axis can't belong to R[i]. Try
#Try:
#GmiRes([5,2,3],2,[{1},{1,2},{3}] );
GmiRes:=proc(m,i,R) local i1,d,j,S,s,k:
option remember:
d:=nops(m):

if nops(R)<>d then
 RETURN(FAIL):
fi:

if not (1<=i and i<=d) then
 RETURN(FAIL):
fi:


if min(seq(m[i1]-m[i1+1],i1=1..d-1))<0 then
 RETURN(0):
fi:


if {seq(m[i1],i1=1..i-1),seq(m[i1],i1=i+1..d)}={0} then
if member(m[i],R[i]) then
 RETURN(0):
else
 RETURN(1):
fi:

fi:

S:={seq(j,j=1..m[i])}  minus R[i]:

add(add(GmiRes([op(1..i-1,m), m[i]-s,op(i+1..d,m)],k,R),s in S),k=1..i-1)+
add(add(GmiRes([op(1..i-1,m), m[i]-s,op(i+1..d,m)],k,R), s in S),k=i+1..d):
end:


#GmRes(m,R): inputs a list m of non-negative integers m, indicating a point in (let d:=nops(m)) the d-dimensional lattice
#and a list of sets R of length d
# outputs the number of positive lattice walks in x[1]>=..x[d]>=0
#from the origin to m such that no runs in the x[i] direction can
#have a length that belongs to R[i]. Try:
#GmRes([4,5,3],[{1},{2},{1}]);
GmRes:=proc(m,R) local i:
option remember:
add(GmiRes(m,i,R),i=1..nops(m)):
end:

#Fseq(RES,N): inputs a list R (let d:=nops(R)) of sets of positive integers outputs the sequence [seq(FmRes([n$d],RES), n=1..N)]
Fseq:=proc(RES,N) local d,n:
d:=nops(RES):
[seq(FmRes([n$d],RES), n=1..N)]:
end:


#Gseq(R,N): inputs a list R (let d:=nops(R)) of sets of positive integers outputs the sequence [seq(GmRes([n$d],R), n=1..N)]
Gseq:=proc(R,N) local d,n:
d:=nops(R):
[seq(GmRes([n$d],R), n=1..N)]:
end:

#YT(L): The set of Young Tableaux of shape L. Try:
#YT([2,2,2]);
YT:=proc(L) local i1,k,gu,n,L1,mu,mu1,i:
option remember:
n:=convert(L, `+`):
k:=nops(L):
if L=[] then
 RETURN({[]}):
fi:

if min(seq(L[i1]-L[i1+1],i1=1..nops(L)-1))<0 then
 RETURN({}):
fi:


gu:={}:

for i from 1 to  k-1 do

 if L[i]>L[i+1] then
  L1:=[op(1..i-1,L),L[i]-1,op(i+1..nops(L),L)]:
  mu:=YT(L1):

 gu:=gu union {seq([op(1..i-1,mu1),[op(mu1[i]),n],op(i+1..k,mu1)],mu1 in mu) } :
fi:
od:

if L[k]>=2 then
  L1:=[op(1..k-1,L),L[k]-1]:
   mu:=YT(L1):
  gu:=gu union {seq([op(1..k-1,mu1),[op(mu1[k]),n]],mu1 in mu) } :
 else
   L1:=[op(1..k-1,L)]:
    mu:=YT(L1):
  gu:=gu union {seq([op(1..k-1,mu1),[n]],mu1 in mu) } :
 fi:
gu:
end:


#MacM(a): MacMahon's formula for the number of Young tableaux of shape a
MacM:=proc(a) local k,i,j:
k:=nops(a):
convert(a,`+`)!/mul((a[i]+k-i)!,i=1..nops(a))*
mul(mul(a[i]-a[j]+j-i,j=i+1..k),i=1..k):
end:


#Runs1(L): Given  an incerasing list of postive integers outputs the sequence of runs. Try
#Runs1([5,7,8,10]);
Runs1:=proc(L) local i:
if L=[] then
 RETURN([]):
fi:

if nops(L)=1 then
 RETURN([1]):
fi:

for i from 1 to nops(L)-1 while L[i+1]-L[i]=1 do od:
[i,op(Runs1([op(i+1..nops(L),L)]))]:
end:


#GmResSet(m,R): The set of Young Tableau of shaper m such that the i-th rown has no run-lenghths in R[i]. Try:
#GmResSet([4,3,2],[{1},{1},{1}]);
GmResSet:=proc(m,R) local j,gu,mu,mu1:

mu:=YT(m):

gu:={}:

for mu1 in mu do
  if {seq(evalb({op(Runs1(mu1[j]))} intersect R[j]={}),j=1..nops(mu1))}={true} then
   gu:=gu union {mu1}:
  fi:
od:

gu:

end:


#Fkx(k,x): The generating function of all words in {1,...,k}. Try:
#Fkx(3,x);
Fkx:=proc(k,x) local F,i,var,eq,i1:

 var:={seq(F[i],i=1..k)}:
 eq:={}:

 for i from 1 to k do

   eq:=eq union {F[i]=x[i]/(1-x[i])*(1+add(F[i1],i1=1..k)-F[i])}:

 od:

var:=solve(eq,var):

normal(1+add(subs(var,F[i]),i=1..k)):

end:


#FkxRes(k,R,x): The generating function of all words in {1,...,m} with no i-runs of length in the set R[i].
#Try:
#FkxRes(k,[{1},{2}],x);
FkxRes:=proc(k,R,x) local F,i,var,eq,i1,lu,j:
option remember:
 var:={seq(F[i],i=1..k)}:
 eq:={}:

 for i from 1 to k do

 lu:=x[i]/(1-x[i])- add(x[i]^j,j in R[i]):
  eq:=eq union {F[i]=lu*(1+add(F[i1],i1=1..k)-F[i])}:
 od:

var:=solve(eq,var):

normal(1+add(subs(var,F[i]),i=1..k)):

end:




#FmResViaGF(m,R): Same as FmRes(m,R) (q.v.) but via Generating functions. For checking purposes only
#and a list of sets R of length d
# outputs the number of positive lattice walks from the origin to m such that no runs in the x[i] direction can
#have a length that belongs to R[i]. Try:
#FmResViaGF([4,5,3],[{1},{2},{1}]);
FmResViaGF:=proc(m,R) local gu,x,i:
gu:=FkxRes(nops(m),R,x):


for i from 1 to nops(m) do
 gu:=normal(coeff(taylor(gu,x[i]=0,m[i]+1),x[i],m[i])):
od:

gu:

end:


###start symbolic in r

#FmiResR(m,i,R,r): inputs a list m of non-negative integers m, indicating a point in (let d:=nops(m)) the d-dimensional lattice
#and 1<=i<=d,  and also a list of length d, of of linear-affine expressions in the variable r,
#outputs the number of lattice walks that end in a run of steps parallel to the i-th axis.
#and that each consecutive run in the i-th axis can't belong to the range of R[i]. Try
#Try:
#FmiResR([5,2,3],2,[2*r+2,2*r+2,2*r+2],r );
FmiResR:=proc(m,i,R,r) local i1,d,j,S,s,k,lu,r1:
option remember:
d:=nops(m):

if nops(R)<>d then
 RETURN(FAIL):
fi:

if not (1<=i and i<=d) then
 RETURN(FAIL):
fi:


if ({seq(m[i1],i1=1..i-1),seq(m[i1],i1=i+1..d)}={0} or {seq(m[i1],i1=1..i-1),seq(m[i1],i1=i+1..d)}={}) then
lu:={seq(subs(r=r1,R[i]),r1=0..m[i])}:
if member(m[i],lu) then
 RETURN(0):
else
 RETURN(1):
fi:

fi:

lu:={seq(subs(r=r1,R[i]),r1=0..m[i])}:
S:={seq(j,j=1..m[i])}  minus lu:

add(add(FmiResR([op(1..i-1,m), m[i]-s,op(i+1..d,m)],k,R,r),s in S),k=1..i-1)+
add(add(FmiResR([op(1..i-1,m), m[i]-s,op(i+1..d,m)],k,R,r), s in S),k=i+1..d):

end:


#FmResR(m,R,r): inputs a list m of non-negative integers m, indicating a point in (let d:=nops(m)) the d-dimensional lattice
#and a list of sets R of length d of sets of linear-affine expressions in the variable r,
# outputs the number of positive lattice walks from the origin to m such that no runs in the x[i] direction can
#have a length that belongs to the range of R[i]. Try:
#FmResR([4,5,3],[2*r+2,2*r+2,2*r+2]);
FmResR:=proc(m,R,r) local i:
option remember:
add(FmiResR(m,i,R,r),i=1..nops(m)):
end:



#FseqR(R,r,N): inputs a list R (let d:=nops(R)) of affine linear expressions in the variable r,
# outputs the sequence [seq(FmResR([n$d],R,r), n=1..N)]. Try:
#FseqR([2*r+2,2*r+2,2*r+2],r,30);
FseqR:=proc(R,r,N) local d,n:
d:=nops(R):
[seq(FmResR([n$d],R,r), n=1..N)]:
end:



#FkxResR(k,R,r,x): The generating function of all words in {1,...,m} with no i-runs of length in the range of R[i] where R[i]
#is an affine-linear expression in r
#Try:
#FkxResR(k,[2*r+1,2*r+1],x);
FkxResR:=proc(k,R,r,x) local F,i,var,eq,i1,lu,j:
option remember:
 var:={seq(F[i],i=1..k)}:
 eq:={}:

 for i from 1 to k do
 lu:=normal(x[i]/(1-x[i])-x[i]^coeff(R[i],r,0)/(1-x[i]^coeff(R[i],r,1))):
  eq:=eq union {F[i]=lu*(1+add(F[i1],i1=1..k)-F[i])}:
 od:

var:=solve(eq,var):

normal(1+add(subs(var,F[i]),i=1..k)):

end:




#FmResRViaGF(m,R,r): Same as FmResR(m,R,r) (q.v.) but via Generating functions. For checking purposes only
#and a list of sets R of length d
# outputs the number of positive lattice walks from the origin to m such that no runs in the x[i] direction can
#have a length that belongs to R[i]. Try:
#FmResRViaGF([4,5,3],[2*r+1,3*r+1,3*r+2]);
FmResRViaGF:=proc(m,R,r) local gu,x,i:
gu:=FkxResR(nops(m),R,r,x):


for i from 1 to nops(m) do
 gu:=normal(coeff(taylor(gu,x[i]=0,m[i]+1),x[i],m[i])):
od:

gu:

end:




#GmiResR(m,i,R,r): inputs a list m of non-negative integers m, indicating a point in (let d:=nops(m)) the d-dimensional lattice
#and 1<=i<=d,  and also a list of length d, of affine-linear expressions in the variable r
#outputs the number of lattice walks  in x[1]>=..x[d]>=0
#that end in a run of steps parallel to the i-th axis.
#and that each consecutive run in the i-th axis can't belong to the range of R[i]. Try

GmiResR:=proc(m,i,R,r) local i1,d,j,S,s,k,lu,r1:
option remember:
d:=nops(m):

if nops(R)<>d then
 RETURN(FAIL):
fi:

if not (1<=i and i<=d) then
 RETURN(FAIL):
fi:


if min(seq(m[i1]-m[i1+1],i1=1..d-1))<0 then
 RETURN(0):
fi:

if ({seq(m[i1],i1=1..i-1),seq(m[i1],i1=i+1..d)}={0} or {seq(m[i1],i1=1..i-1),seq(m[i1],i1=i+1..d)}={}) then
lu:={seq(subs(r=r1,R[i]),r1=0..m[i])}:
if member(m[i],lu) then
 RETURN(0):
else
 RETURN(1):
fi:

fi:

lu:={seq(subs(r=r1,R[i]),r1=0..m[i])}:
S:={seq(j,j=1..m[i])}  minus lu:

add(add(GmiResR([op(1..i-1,m), m[i]-s,op(i+1..d,m)],k,R,r),s in S),k=1..i-1)+
add(add(GmiResR([op(1..i-1,m), m[i]-s,op(i+1..d,m)],k,R,r), s in S),k=i+1..d):

end:


#GmResR(m,R,r): inputs a list m of non-negative integers m, indicating a point in (let d:=nops(m)) the d-dimensional lattice
#and a list of linear affine expressions of r  of length d
# outputs the number of positive lattice walks in x[1]>=..x[d]>=0
#from the origin to m such that no runs in the x[i] direction can
#have a length that belongs to the range of R[i]. Try:
#GmResR([4,5,3],[2*r+1,2*r+1,3*r+1],r);
GmResR:=proc(m,R,r) local i:
option remember:
add(GmiResR(m,i,R,r),i=1..nops(m)):
end:


#GseqR(R,r,N): inputs a list R (let d:=nops(R)) of affine linear expressions in the variable r,
# outputs the sequence [seq(GmResR([n$d],R,r), n=1..N)]. Try:
#GseqR([2*r+2,2*r+2,2*r+2],r,30);
GseqR:=proc(R,r,N) local d,n,i:

if not type(R,list) then
 print(`Bad input`):
 RETURN(FAIL):
fi:

if not {seq(degree(R[i],r),i=1..nops(R))} minus {0,1}={} then
 print(`Bad input`):
 RETURN(FAIL):
fi:

if not {seq(type(coeff(R[i],r,1),integer),i=1..nops(R))}={true} then
 print(`Bad input`):
 RETURN(FAIL):
fi:

if not {seq(evalb(coeff(R[i],r,1)>=0),i=1..nops(R))}={true} then
 print(`Bad input`):
 RETURN(FAIL):
fi:

if not {seq(evalb(coeff(R[i],r,0)>=0),i=1..nops(R))}={true} then
 print(`Bad input`):
 RETURN(FAIL):
fi:


d:=nops(R):
[seq(GmResR([n$d],R,r), n=1..N)]:
end:



#InfoF(RES,N0,n,a,MaxC,N1): Inputs
#(i) a list RES of sets of posivive integers, let k=nops(RES). outputs
#(ii) a fairly large pos. integer N0
#(iii) symbols n and a
#(iv) a pos. integer MaxC
#(v): A large positive integer N1
#Outputs a list consisting of
#(i) The first N0 terms of the sequence "number of walks in k-dim hyper-cubic lattice
#from the origin to [n,n,n] without runs in the i-th direction of length in R[i]
#(ii) A linear recurrence equation with poly coefficients of complexity <=MaxC
#satisfied by the sequence in terms of a(n), followed by the initial conditions . If none found it only returns the above list
#(iii) The value of the N1-th term of the sequence. Try
#InfoF([{1},{1}],50,n,a,10,100);

InfoF:=proc(RES,N0,n,a,MaxC,N1) local N,gu,ope,i1,lu:

gu:=Fseq(RES,N0):

ope:=Findrec(gu,n,N,MaxC):

if ope=FAIL then
 RETURN([gu]):
fi:

lu:=SeqFromRec(ope,n,N,[op(1..degree(ope,N),gu)],N1)[N1]:

[gu,[add(coeff(ope,N,i1)*a(n+i1),i1=0..degree(ope,N))=0,{seq(a[i1]=gu[i1],i1=1..degree(ope,N))}],lu]:

end:



#InfoFr(RES,r,N0,n,a,MaxC,N1): Inputs
#(i) a list RES of sets of affine linear expressions in r, let k=nops(RES). outputs
#(ii) a fairly large pos. integer N0
#(iii) symbols n and a
#(iv) a pos. integer MaxC
#(v): A large positive integer N1
#Outputs a list consisting of
#(i) The first N0 terms of the sequence "number of walks in k-dim hyper-cubic lattice
#from the origin to [n,n,n] without runs in the i-th direction of length in  the range of R[i]
#(ii) A linear recurrence equation with poly coefficients of complexity <=MaxC
#satisfied by the sequence in terms of a(n), followed by the initial conditions . If none found it only returns the above list
#(iii) The value of the N1-th term of the sequence. Try
#InfoFr([r+3,r+3],r,50,n,a,10,100);

InfoFr:=proc(RES,r,N0,n,a,MaxC,N1) local N,gu,ope,i1,lu:

gu:=FseqR(RES,r,N0):

ope:=Findrec(gu,n,N,MaxC):

if ope=FAIL then
 RETURN([gu]):
fi:

lu:=SeqFromRec(ope,n,N,[op(1..degree(ope,N),gu)],N1)[N1]:

[gu,[add(coeff(ope,N,i1)*a(n+i1),i1=0..degree(ope,N))=0,{seq(a[i1]=gu[i1],i1=1..degree(ope,N))}],lu,ope]:

end:




#InfoG(RES,N0,n,a,MaxC,N1): Inputs
#(i) a list RES of sets of posivive integers, let k=nops(RES). outputs
#(ii) a fairly large pos. integer N0
#(iii) symbols n and a
#(iv) a pos. integer MaxC
#(v): A large positive integer N1
#Outputs a list consisting of
#(i) The first N0 terms of the sequence "number of walks in k-dim hyper-cubic lattice
#from the origin to [n$k] staying in x[1]>= >=x[k] without runs in the i-th direction of length in R[i]
#(ii) A linear recurrence equation with poly coefficients of complexity <=MaxC
#satisfied by the sequence in terms of a(n), followed by the initial conditions . If none found it only returns the above list
#(iii) The value of the N1-th term of the sequence. Try
#InfoG([{1},{1}],50,n,a,10,100);

InfoG:=proc(RES,N0,n,a,MaxC,N1) local N,gu,ope,i1,lu:

gu:=Gseq(RES,N0):

ope:=Findrec(gu,n,N,MaxC):

if ope=FAIL then
 RETURN([gu]):
fi:

lu:=SeqFromRec(ope,n,N,[op(1..degree(ope,N),gu)],N1)[N1]:

[gu,[add(coeff(ope,N,i1)*a(n+i1),i1=0..degree(ope,N))=0,{seq(a[i1]=gu[i1],i1=1..degree(ope,N))}],lu]:

end:



#InfoGr(RES,r,N0,n,a,MaxC,N1): Inputs
#(i) a list RES of sets of affine linear expressions in r, let k=nops(RES). outputs
#(ii) a fairly large pos. integer N0
#(iii) symbols n and a
#(iv) a pos. integer MaxC
#(v): A large positive integer N1
#Outputs a list consisting of
#(i) The first N0 terms of the sequence "number of walks in k-dim hyper-cubic lattice
#from the origin to [n,..,n] staying in x[1]>=x[2]>=... x[k] without runs in the i-th direction of length in  the range of R[i]
#(ii) A linear recurrence equation with poly coefficients of complexity <=MaxC
#satisfied by the sequence in terms of a(n), followed by the initial conditions . If none found it only returns the above list
#(iii) The value of the N1-th term of the sequence. Try
#InfoGr([r+3,r+3],r,50,n,a,10,100);

InfoGr:=proc(RES,r,N0,n,a,MaxC,N1) local N,gu,ope,i1,lu:

gu:=GseqR(RES,r,N0):

ope:=Findrec(gu,n,N,MaxC):

if ope=FAIL then
 RETURN([gu]):
fi:

lu:=SeqFromRec(ope,n,N,[op(1..degree(ope,N),gu)],N1)[N1]:

[gu,[add(coeff(ope,N,i1)*a(n+i1),i1=0..degree(ope,N))=0,{seq(a[i1]=gu[i1],i1=1..degree(ope,N))}],lu,ope]:

end:


#ZinnPlus(resh,n): An estimate for the asympotics of resh[n] in the form C*mu^n*n^theta. Try:
#L:=Gseq([{1},{1}],100): ZinnPlus(L,n);
ZinnPlus:=proc(resh,n) local gu,C,lu,i,err,d:
gu:=Zinn(resh):

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

lu:=n^gu[1]*gu[2]^n:

C:=[seq(evalf(resh[i]/subs(n=i,lu)),i=nops(resh)-4..nops(resh))]:

err:=abs(C[1]-C[nops(C)]):

if err>1/100 then
 RETURN(FAIL):
fi:

d:=trunc(-log[10](err)):

C:=evalf(C[nops(C)],d):
C*lu

end:


#MyAsy(L,n,k): tries to fit the sequence L to be of the form  C*mu^n*n^theta*(1+c1/n+...+ck/n^k). Try:
#L:=Gseq([{1},{1}],50); MyAsy(L,n,3);

MyAsy:=proc(L,n,k) local gu,logC,logmu,theta,c,i,eq,var,C,mu,x,ku:

if nops(L)<50 then
 print(`List must be at least of length 50`):
 RETURN(FAIL):
fi:

if 10+k>nops(L) then
 print(`k must be at most`, nops(L)-10 ):
fi:

gu:=logC+logmu*n+theta*log(n)+add(c[i]/n^i,i=1..k):

var:={logC,logmu,theta,seq(c[i],i=1..k)}:

eq:={seq(log(L[i])-subs(n=i,gu),i=nops(L)-k-2..nops(L))}:

var:=evalf(solve(eq,var)):

C:=exp(subs(var,logC)):
mu:=exp(subs(var,logmu)):
theta:=subs(var,theta):

C:=evalf(C,10):
mu:=evalf(mu,10):
theta:=evalf(theta,10):

ku:= exp(subs(var,add(c[i]*x^i,i=1..k))):
ku:=taylor(ku,x=0,k+1):
ku:=add(evalf(coeff(ku,x,i),10)/n^i,i=0..k):
evalf(C*mu^n*n^theta*ku,10):
end:





#GseqVA(R,N,n): Verbose form of Gseq(R,N) (q.v.) and also gives the asympototics to order 2. N has to be at least 50
#Try:
#GseqVA([{1},{1},{1}],50,n);
GseqVA:=proc(R,N,n) local L,d,t0,i:

t0:=time():
if not type(R,list) then
 print(`Bad input`):
 RETURN(FAIL):
fi:

if not type(N,integer) and N>=50 then
 print(N, `must be an integer >=50`):
 RETURN(FAIL):
fi:


L:=Gseq(R,N):

d:=nops(R):

print(`The sequence enumerating Young tableaux of shape`, [n$d], `such that `):

for i from 1 to nops(R) do
 print(`The `, i, ` -th row does not have runs in the set `, R[i]):
 print(``):
od:

print(`for n from 1 to` , N,  ` is `):
print(``):
print(L):
print(``):
print(`The estimated asymptotics is`):
print(``):
print(MyAsy(L,n,2)):
print(``):

print(`------------------`):
print(``):
print(`This took `, time()-t0, `seconds. `):
print(``):
end:



#GseqRVA(R,r,N,n): Verbose form of GseqR(R,r,N) (q.v.) and also gives the asympototics to order 2. N has to be at least 50
#Try:
#GseqRVA([2*r,2*r,2*r],r,50,n);
GseqRVA:=proc(R,r,N,n) local L,d,t0,i:

t0:=time():
if not type(R,list) then
 print(`Bad input`):
 RETURN(FAIL):
fi:

if not (type(N,integer) and N>=50 )then
 print(N, `must be an integer >=50`):
 RETURN(FAIL):
fi:


L:=GseqR(R,r,N):

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

d:=nops(R):

print(`The sequence enumerating Young tableaux of shape`, [n$d], `such that `):

for i from 1 to nops(R) do
 print(`The `, i, ` -th row does not have runs in the arithmetical progression (starting at r=0) `, R[i]):
 print(``):
od:

print(`for n from 1 to` , N,  ` is `):
print(``):
print(L):
print(``):
print(`The estimated asymptotics is`):
print(``):
print(MyAsy(L,n,2)):
print(``):

print(`------------------`):
print(``):
print(`This took `, time()-t0, `seconds. `):
print(``):
end:



#MyAsyC(L,n,k,theta,mu): tries to fit the sequence L to be of the form  C*mu^n*n^theta*(1+c1/n+...+ck/n^k). 
#when mu and theta are given. Try:
#L:=Gseq([{},{}],50); MyAsyC(L,n,3,-3/2,4);
MyAsyC:=proc(L,n,k,theta,mu) local gu,logC,logmu,c,i,eq,var,C,x,ku:

if nops(L)<50 then
 print(`List must be at least of length 50`):
 RETURN(FAIL):
fi:

if 10+k>nops(L) then
 print(`k must be at most`, nops(L)-10 ):
fi:

logmu:=log(mu):
gu:=logC+logmu*n+theta*log(n)+add(c[i]/n^i,i=1..k):

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

eq:={seq(log(L[i])-subs(n=i,gu),i=nops(L)-k..nops(L))}:

var:=evalf(solve(eq,var)):

C:=exp(subs(var,logC)):
C:=evalf(C,10):
ku:= exp(subs(var,add(c[i]*x^i,i=1..k))):
ku:=taylor(ku,x=0,k+1):
ku:=add(evalf(coeff(ku,x,i),10)/n^i,i=0..k):
evalf(C*mu^n*n^theta*ku,10):
end: