Help:=proc(): print( ` GuessH(L,n,N) , SpellOut(ope,n,N,f)`): print(`HtoSeq(ope,n,N,Ini,M)`): end: #old stuff from C13.txt Help13:=proc(): print(` LIF(P,u,N), GuessH1(L,ORD,DEG,n,N)`): end: #Lagrange Inversion applied to tree-enumeration #of ordered #Starting the holonomic ansatz (a.k.a P-recursive) #LIF: Inputs an expression P of the variable #u that possesses #a Maclaurin expansion in u, and outputs the #first N terms of the series expansion of the #unique f.p.s u(t) satisfying the functional eq. #u(t)=t*P(u(t)). For example, #LIF(1+u^2,u,10); should return something with #Catalan numbers LIF:=proc(P,u,N) local n: [seq(coeff(taylor(P^n, u=0,n),u,n-1)/n,n=1..N)]: end: #GuessH1(L,ORD,DEG,n,N): Inputs a sequence of numbers L #and pos. integer ORD and non-neg. integer DEG and letters #n and N and outputs an linear recurrence operator with #poly. coeffs. (conj.) annihilating the sequence L #nops(L)>=(1+ORD)*(1+DEG)+ORD+6 GuessH1:=proc(L,ORD,DEG,n,N) local ope,a,i,j,var,eq,n1,var1: if nops(L)<(1+ORD)*(1+DEG)+ORD+6 then print(`We need at least`, (1+ORD)*(1+DEG)+ORD+6, `terms in L`): RETURN(FAIL): fi: ope:=add(add(a[i,j]*n^i*N^j,i=0..DEG),j=0..ORD): var:={seq(seq(a[i,j],i=0..DEG),j=0..ORD)}: eq:={seq(add(subs(n=n1, coeff(ope,N,j))*L[n1+j],j=0..ORD), n1=1..(1+DEG)*(1+ORD)+6)}: var1:=solve(eq,var): ope:=subs(var1,ope): if ope=0 then RETURN(FAIL): fi: ope:=subs({seq(v=1, v in var)}, ope): ope:=add(factor(coeff(ope,N,j))*N^j, j=0..degree(ope,N)): eq:=expand({seq(add(subs(n=n1, coeff(ope,N,j))*L[n1+j],j=0..ORD), n1=1..nops(L)-ORD)}): if eq<>{0} then print(ope, `didn't work out`): RETURN(FAIL): else RETURN(ope): fi: end: #end old stuff #start new stuff #GuessH(L,n,N): Inputs a sequence of numbers L # and symbols n and N #and finds a minimal (in the sense of ord+deg) # linear recurrence operator with #poly. coeffs. (conj.) annihilating the sequence L with #nops(L)>=(1+ORD)*(1+DEG)+ORD+6 GuessH:=proc(L,n,N) local K,ope,DEG,ORD: ORD:=1: DEG:=0: for K from 1 while nops(L)>=(1+ORD)*(1+DEG)+ORD+6 do for ORD from 1 to K do DEG:=K-ORD: ope:=GuessH1(L,ORD,DEG,n,N): if ope<>FAIL then RETURN(ope): fi: od: od: FAIL: end: #SpellOut(ope,n,N,f) #writes the fact ope(n,N)f(n)=0 in common languae SpellOut:=proc(ope,n,N,f) local i: add(coeff(ope,N,i)*f(n+i),i=0..degree(ope,N))=0: end: #HtoSeq(ope,n,N,Ini,M): Inputs a linear recurrence operator #ope phrased in terms of the discrete variable n and #its corresponding shift operator N, and ORD=degree(ope,N) #values for the sequence outputs the first M terms of #the sequence. For example, #HtoSeq(n+1-N,n,N,[1],4); should yield #[1,2,6,24] #N^(-ORD)*ope(N,n) HtoSeq:=proc(ope,n,N,Ini,M) local L,n1,ORD,BenMiles,kirsten,i: ORD:=degree(ope,N): if nops(Ini)<>ORD then RETURN(FAIL): fi: BenMiles:=expand(subs(n=n-ORD,ope)/N^ORD): L:=Ini: for n1 from nops(Ini)+1 to M do kirsten:=subs(n=n1,coeff(BenMiles,N,0)): if kirsten=0 then print(`blows up at`, n1): print(`so far we have`, L): RETURN(FAIL): fi: L:=[op(L), -add(subs(n=n1,coeff(BenMiles,N,-i))*L[n1-i],i=1..ORD)/kirsten]: od: L: end: