## Kristen Lew, Homework 5 ## Exp Math Spring 2012 # Use whatever you like Help:=proc(): print(` Old processes: Di(L), GP(L,x), GQP1(L,x,m), GQP(L,x)`): print(` New processes: GFtoSeq(f,q,m), pmn(m,n), F(m,x) `): end: ### Problem 2 ## GFtoSeq(f,q,m): inputs an expression f in the variable q # and outputs the list L of number of length m such that L[i] is the coefficient # of q^i in the Maclaurin expansion of f with respect to variable q for i from 1 to m. GFtoSeq:=proc(f,q,m) local T,L,i: T:=taylor(f,q=0,m+1): L:=[]: for i from 1 to m do L:=[op(L), coeff(T,q,i)]: od: L: end: ### Problem 3 ## pmn(m,n): inputs a positive intger m and outputs a quasipolynomial of degree m # and order lcm(1,…,m) that describes p_m(n). Test program for m=1,2,3,4,5. pmn:=proc(m,n) local f, lcm1, F, i: f:=1/(product(1-(q^i), i=1..m)): lcm1:=lcm(seq(i,i=1..m)): F:=GFtoSeq(f,q,lcm1*(10)): GQP1(F,n,lcm1): end: ## GQP(L,x) inputs a list of number L and a variable x 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: ## GQP1(L,x,m) inputs a list of numbers L, a variable x, # and a positive integer m, and outputs a quasi-polynomial, Q, # of order m describing a fitting L such that if # n is i mod m then, L[n] equals the Q[i](n) GQP1:=proc(L,x,m) local Q,i,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: # 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: Di:=proc(L) local i: [ seq( L[i+1]-L[i],i=1..nops(L)-1)]: end: ### Problem 4 ## F(m,x) that expresses trunc(x/m) as a quasi-polynomial in x of order m for integer x. F:=proc(m,x) local T: T:= [seq(trunc(i/m),i=1..m*(10))]: GQP1(T,x,m): end: ### Problem 5 ## > T:= n -> round(n^2/12) - trunc(n/4)*trunc((n+2)/4); # T := proc (n) options operator, arrow; round((1/12)*n^2)- # trunc((1/4)*n)*trunc((1/4)*n+1/2) end proc ## > Tn:=[seq(T(i), i=0..100)]; # Tn := [0, 0, 0, 1, 0, 1, 1, 2, 1, 3, 2, 4, 3, 5, 4, 7, 5, 8, 7, 10, 8, 12, 10, # 14, 12, 16, 14, 19, 16, 21, 19, 24, 21, 27, 24, 30, 27, 33, 30, 37, 33, 40, 37, # 44, 40, 48, 44, 52, 48, 56, 52, 61, 56, 65, 61, 70, 65, 75, 70, 80, 75, 85, 80, # 91, 85, 96, 91, 102, 96, 108, 102, 114, 108, 120, 114, 127, 120, 133, 127, 140, # 133, 147, 140, 154, 147, 161, 154, 169, 161, 176, 169, 184, 176, 192, 184, 200, # 192, 208, 200, 217, 208] ## GQP(Tn,x); # [(1/48)*x^2-(1/24)*x+1/48, (1/48)*x^2+(1/12)*x-1/4, (1/48)*x^2-(1/24)*x-1/16, # 1/3+(1/12)*x+(1/48)*x^2, (1/48)*x^2-(1/24)*x-5/16, (1/48)*x^2+(1/12)*x-1/4, # 13/48-(1/24)*x+(1/48)*x^2, (1/48)*x^2+(1/12)*x, (1/48)*x^2-(1/24)*x-5/16, # (1/48)*x^2+(1/12)*x+1/12, (1/48)*x^2-(1/24)*x-1/16, (1/48)*x^2+(1/12)*x]