Help:=proc(): print(` GP(L,x),GQP1(L,x,m), GQP(L,x) `): end: Di:=proc(L) local i: [ seq( L[i+1]-L[i],i=1..nops(L)-1)]: end: #GP(L,x): The unique polynomial P(x) #expressed in terms of binomial(x,k) #such that P(i)=L[i+1] for i from 0 to nops(L)-1 GP:=proc(L,x) local i,L1,P: L1:=L: P:=0: for i from 1 to nops(L)-2 do P:=P+L1[1]*binomial(x,i-1): L1:=Di(L1): if convert(L1,set)={0} then RETURN(factor(expand( subs(x=x-1,P)))): fi: od: FAIL: end: #GQP1(L,x,m): inputs a list of numbers L, a variable x #and A pos. integer m, and outputs a quasi-poly of order m, #let's call it Q #fitting L such that if n is i mod m then L[n] equals #Q[i](n) GQP1:=proc(L,x,m) local Q,i,j,L1,P1: Q:=[]: for i from 1 to m do L1:=[ seq(L[i+m*j],j=0..trunc((nops(L)-i)/m) ) ]: P1:=GP(L1,x): if P1=FAIL then RETURN(FAIL): else P1:=expand(subs(x=(x-i)/m+1,P1)): Q:=[ op(Q), P1]: fi: od: Q: end: #GQP(L,x): inputs a list of numbers and tries to guess # a quasi-polynomial fitting the data, or returns FAIL GQP:=proc(L,x) local m, Ben: for m from 1 to nops(L)/2 do Ben:=GQP1(L,x,m): if Ben<>FAIL then RETURN(Ben): fi: od: FAIL: end: