#C8.txt: Feb. 18, 2019 Expmath, Dr. Z. Help:=proc(): print(` CtoSeq(C,N), GuessPolC(C,N,x,r)`): end: ##stuff from C7.txt Help7:=proc(): print(` ProbConts(), Wt(L,x) , GFp(LL,x) , Comp(n,k) , CompU(n,k) `): print(`GuessMulPol1(L,x,k,d), GuessMulPol(L,x,k) `): end: #from C4.txt #Comp(n,k): inputs non-neg. integers n and k and outputs the SET of # compositions of n into k parts (vectors of non-neg. integers that add-up to n) #of length k Comp:=proc(n,k) local S,a1,S1, s1: option remember: if k<0 or n<0 then RETURN({}): fi: if n=0 then RETURN({[0\$k]}): fi: if k=0 then if n=0 then RETURN({[]}): else RETURN({}): fi: fi: S:={}: for a1 from 0 to n do S1:=Comp(n-a1,k-1): S:=S union {seq([a1,op(s1)], s1 in S1)}: od: S: end: #end C4.txt #CompU(n,k) the set of vectors with k non-neg. componets that add up to <=n CompU:=proc(n,k) local n1: {seq(op(Comp(n1,k)),n1=0..n)}: end: #GuessMulPol1(L,x,k,d): inputs a list L of data of the form [pt,value] where #pt is a list of length k, and value is the value. Guesses a polynomial #in x[1], ..., x[k] of degree d,P, such that P(L[i][1])=L[i][2] GuessMulPol1:=proc(L,x,k,d) local S,a,s,P,i,var,eq,i1: if not ( type(L,list) and {seq(type(L[i],list),i=1..nops(L))}={true} and {seq(nops(L[i]),i=1..nops(L))}={2} and {seq(type(L[i][1],list),i=1..nops(L))}={true} and {seq(nops(L[i][1]),i=1..nops(L))}={k}) then print(`bad input`): RETURN(FAIL): fi: S:=CompU(d,k): P:=add(a[s]*mul(x[i]^s[i],i=1..k), s in S): var:={seq(a[s],s in S)}: eq:={seq(subs({seq(x[i1]=L[i][1][i1],i1=1..k)},P)=L[i][2],i=1..nops(L))}: if nops(eq)-nops(var)<=3 then print(`not enough data `): RETURN(FAIL): fi: var:=solve(eq,var): if var=NULL then RETURN(FAIL): fi: P:= subs(var,P): if P=0 then RETURN(FAIL): else RETURN(P): fi: end: #Added after class ended #GuessMulPol(L,x,k): inputs a list L of data of the form [pt,value] where #pt is a list of length k, and value is the value. Guesses a polynomial #in x[1], ..., x[k] of degree <=Maxd ,P, such that P(L[i][1])=L[i][2] GuessMulPol:=proc(L,x,k) local d,P,d1: for d from 1 while nops(L)-nops(CompU(d,k))>3 do od: d:=d-1: for d1 from 1 to d do P:=GuessMulPol1(L,x,k,d1): if P<>FAIL then RETURN(P): fi: od: FAIL: end: ### #CtoSeq(C,N): inputs a C-finite sequence C=[INI,Rec] and a pos. integer N and #outputs the first N terms starting at n=0. CtoSeq:=proc(C,N) local L,INI,R,nt,i: INI:=C[1]: R:=C[2]: if not (type(C,list) and nops(C)=2 and type(C[1],list) and type(C[2],list) and nops(C[1])=nops(C[2])) then print(`bad input`): RETURN(FAIL): fi: L:=INI: while nops(L)<=N do nt:=add(R[i]*L[-i],i=1..nops(R)): L:=[op(L),nt]: od: L: end: #GuessPolC(C,N,x,r): Inputs a C-finite sequence C=[INI,Rec] (denoting a sequence of order k, say) #k=nops(INI) and a(n)=Rec[1]*a[n-1]+ ... + Rec[k]*a[n-k] . Fib is [[0,1],[1,1]] #and outputs a polynomial P, in x[1], ..., x[r], such that P(a(i),...,a(i+r-1))==0 #N is the number of data gathered GuessPolC:=proc(C,N,x,r) local L1,L2,i,j,P,P0,c: L1:=CtoSeq(C,N): #GuessMulPol(L,x,k) L2:= [ seq( [ [seq(L1[i+j],j=0..r-1)] ,0],i=1..N-r+1)]; P:=GuessMulPol(L2,x,r): if P=FAIL then RETURN(FAIL): fi: P0:=subs({seq(x[i]=i^2,i=1..r)},P): P:=normal(P/P0): P/lcoeff(P): end: