#John Kim
#Use whatever you like


Help:=proc() : print( `IsConst(L),Di(L),IsPol(L) , GuessPol(L,x) `): end:


#Di(L): the sequence of differences of consectutive entries of L
#For example,Di([1,4,9]); should be [3,5].
Di:=proc(L) local i:

[seq(L[i]-L[i-1],i=2..nops(L))]:

end:

#IsConst(L): is the sequence of numbers L a constant sequence:
IsConst:=proc(L) 
   evalb( nops(convert(L,set))<=1 ):
end:


#IsPol(L): inputs a list of length L and returns
#a conjectured degree if it (conj.) a polynomial
#and FAIL otherwise
IsPol:=proc(L) local L1, k:

L1:=L:

for k from 1 to nops(L)-4 do
 L1:=Di(L1):

 if L1=[0$nops(L1)] then
    RETURN(k-1):
 fi:

od:

FAIL:

end:


#GuessPol(L,x): guesses a polynomial expression P(x)
#such that P(i)=L[i] for all i between 1 and nops(L)
GuessPol:=proc(L,x) local d,P,var, eq,a,i,c,Q,start,finish:

start:=time():
d:=IsPol(L):


if d=FAIL then
 finish:=time():
 print(finish-start):
 RETURN(d):
fi:

P:=add( a[i]*x^i, i=0..d):

var:={ seq( a[i], i=0..d)}:

eq:={ seq( subs(x=c,P)=L[c], c=1..d+1)}:


var:=solve(eq,var):


if var=NULL then
  finish:=time():
  print(finish-start):
  RETURN(FAIL):
fi:

Q:=subs( var, P):
if convert([seq(subs(x=c,Q)-L[c], c=d+2..nops(L))],set)={0} then
  finish:=time():
  print(finish-start):
  Q:
else
  finish:=time():
  print(finish-start):
  RETURN(FAIL):
fi:

end:


#GuessRat1(L,x,d1,d2): guesses a rational poly expression P(x)/Q(x)
#such that deg(P)=d1, deg(Q)=d2, and P(i)/Q(i)=L[i] for all i between
#1 and nops(L).  L should be longer than d1+d2+4.

GuessRat1:=proc(L,x,d1,d2) local a,b,P,Q,var,eq,i,c:

if d1+d2+4>nops(L) then
  RETURN(FAIL):
fi:

P:=add(a[i]*x^i, i=0..d1-1)+x^d1:
Q:=add(b[i]*x^i, i=0..d2):

var:={seq(a[i],i=0..d1),seq(b[i],i=0..d2)}:
eq:={seq(subs(x=c,P/Q)=L[c],c=1..nops(L))}:

var:=solve(eq,var):

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

subs(var,P/Q):

end:

#GuessRat(L,x): guesses a rational poly expression P(x)/Q(x)
#that minimizes deg(P)+deg(Q).  Returns FAIL if none is found.

GuessRat:=proc(L,x) local d,d1,d2:

for d from 0 to nops(L)-4 do
  for d1 from 0 to d do
    d2:=d-d1:
    if GuessRat1(L,x,d1,d2) <> FAIL then
      RETURN(GuessRat1(L,x,d1,d2)):
    fi:
  od:
od:

RETURN(FAIL):

end:


#FF(x,k): returns falling factorial of degree k in x.
FF:=proc(x,k) local i:
if k=0 then
  1:
else
  mul(x-i,i=0..k-1):
fi:
end:


#BetterGuessPol(L,x): guesses a polynomial expression P(x)
#such that P(i)=L[i] for all i between 1 and nops(L)
#and FAIL otherwise
BetterGuessPol:=proc(L,x) local L1,k,a,start,finish:

start:=time():
L1:=L:
a[0]:=L1[1]:

for k from 1 to nops(L)-4 do
 L1:=Di(L1):
 a[k]:=L1[1]/k!:

 if L1=[0$nops(L1)] then
    finish:=time():
    print(finish-start):
    RETURN(expand(add(a[i]*FF(x-1,i),i=0..k))):
 fi:

od:

finish:=time():
print(finish-start):
FAIL:

end:

#For the polynomial sequence x^80+x^23+3*i (100 terms), BetterGuessPol runs in 0.032 seconds.
#GuessPol runs in 0.375 seconds.  In this case, over 10 times better!