# Homework 23. 15-April-2012. Pat Devlin # # I don't care about keeping my work private # Help:=proc(): print(`PT(), Mat(), findConway(), RLT1(S,N), RLT2(S,N)`): end proc: HelpAll:=proc(): print(`split(w1, w2), SplitUp(L)`): print(`JC(w), YN(INI,n), JC1(w) `): end proc: ############# # Problem 1 # ############# PT:=proc() local atoms, checkedAtoms, newAtoms, newCheckedAtoms, L, i: atoms:={[3]}: checkedAtoms:={}: while(not (atoms = checkedAtoms)) do newAtoms:={}: newCheckedAtoms:={}: L:=[op(atoms minus checkedAtoms)]: for i from 1 to nops(L) do newAtoms:=newAtoms union {op(SplitUp(JC(L[i])))}: newCheckedAtoms:=newCheckedAtoms union {L[i]}: od: atoms:=atoms union newAtoms: checkedAtoms:=checkedAtoms union newCheckedAtoms: od: return atoms: end proc: ############# # Problem 2 # ############# Mat:=proc() local L, A, i, j: L:=[op(PT())]: A:=[]: for i from 1 to nops(L) do for j from 1 to nops(L) do A:=[op(A), entry(L[i],L[j])]: od: od: return Matrix(nops(L), nops(L), A): end proc: entry:=proc(w1,w2) local L, i, j: L:=SplitUp(JC(w1)): i:=0: for j from 1 to nops(L) do if(L[j] = w2) then i:=i+1: fi: od: return i: end proc: ############# # Problem 3 # ############# with(LinearAlgebra): findConway:=proc() local x, A, p: A:=Mat(): p:=CharacteristicPolynomial(A,x): return FindLargestRoot(p,x): end proc: FindLargestRoot:=proc(p,x) local biggestModulus, biggestIndex, i, L: L:=[fsolve(p=0, x)]: biggestIndex:=1: biggestModulus:=evalf(abs(L[1])): for i from 1 to nops(L) do if (evalf(abs(L[i])) > biggestModulus) then biggestModulus:=evalf(abs(L[i])): biggestIndex:=i: fi: od: return L[biggestIndex]: end proc: ############# # Problem 4 # ############# RLT1:=proc(S,N) local f, u, s: f:=1+ add(u^s / s!, s in (S minus {0})): return LIF(f, u, N): end proc: # RLT1({1},N), RLT1({1,2}, N), and RLT1({1,2,3}, N) are all in Sloane # RLT1({1,3}, N) is not in Sloane, but personally I don't find it particularly interesting ############# # Problem 5 # ############# RLT2:=proc(S,N) local f, u, s: f:=exp(u) - add(u^s / s!, s in (S minus {0})): return LIF(f, u, N): end proc: # RLT2({1}, N) was in Sloane, but not for this tree interpretation # RLT2({1,2}, N), RLT2({1,2,3}, N), and RLT2({1,3}, N) are not in Sloane # (Again, I personally don't find them interesting enough to submit them) ##################### # Stuff from Before # ##################### #################### # From Homework 22 # #################### SplitUp:=proc(L) option remember: local i: for i from 1 to (nops(L)-1) do if(split(L[1..i], L[(i+1)..nops(L)])) then return [op(SplitUp(L[1..i])), op(SplitUp(L[(i+1)..nops(L)]))]: fi: od: return [L]: end proc: ################# # From Class 23 # ################# LIF:=proc(P,u,N) local n: return [seq(coeff(taylor(P^n, u=0,n),u,n-1)/n,n=1..N)]: end proc: ################# # From Class 22 # ################# #The Conway Look-And-Say sequence #JC(w): inputs any list in and outputs its #Yusra transform #a1^c1,a2^c2,..., ar^cr ->[a1,c1,a2,c2,a3,c3, ...] JC:=proc(w) local i,a1,c1,w1,J1: option remember: if w=[] then RETURN([]): fi: a1:=w[1]: for i from 1 to nops(w) while w[i]=a1 do od: c1:=i-1: w1:=w[c1+1..nops(w)]: J1:=JC(w1): [c1,a1,op(J1)]: end: #YN(INI,n): the first n terms of the #sequence : "length of JC^n(INI) YN:=proc(INI,n) local s,i,Cu: s:=[]: Cu:=INI: for i from 1 to n do s:=[op(s), nops(Cu)]: Cu:=JC(Cu): od: s: end: #JC1(w): inputs any list in and outputs its #Yusra transform, not using recursion #a1^c1,a2^c2,..., ar^cr ->[a1,c1,a2,c2,a3,c3, ...] JC1:=proc(w) local i,a1,c1,w1,J1,L: option remember: if w=[] then RETURN([]): fi: L:=[]: w1:=w: while w1<>[] do a1:=w1[1]: for i from 1 to nops(w1) while w1[i]=a1 do od: c1:=i-1: L:=[op(L),c1,a1]: w1:=w1[c1+1..nops(w1)]: od: L: end: #YN1(INI,n): the first n terms of the #sequence : "length of JC^n(INI) YN1:=proc(INI,n) local s,i,Cu: s:=[]: Cu:=INI: for i from 1 to n do s:=[op(s), nops(Cu)]: Cu:=JC1(Cu): od: s: end: ###### #################################### # split(w1,w2) as from the Webpage # #################################### split:=proc(L,R) local n,m: if nops(L)=0 or nops(R)=0 then RETURN(1): fi: n:=op(nops(L),L): m:=op(1,R): if n=m then RETURN(false): fi: if n>=4 and m<=3 then RETURN(true): fi: if n=2 then if m=1 and nops(R)>2 and op(2,R)<>1 and op(3,R)<>op(2,R) then RETURN(true): fi: if m=1 and nops(R)=2 and op(2,R)<>1 then RETURN(true): fi: if m=1 and nops(R)>=3 and op(2,R)=1 and op(3,R)=1 then RETURN(true): fi: if m=3 and nops(R)>=2 and op(2,R)<>3 and not(nops(R)>=4 and op(2,R)=op(3,R) and op(3,R)=op(4,R)) then RETURN(true): fi: if m=3 and nops(R)=1 then RETURN(true): fi: if m>=4 then RETURN(true): fi: fi: if n<>2 then if nops(R)>=5 and op(1,R)=2 and op(2,R)=2 and op(3,R)=1 and op(4,R)<>1 and op(5,R)<>op(4,R) then RETURN(true): fi: if nops(R)=4 and op(1,R)=2 and op(2,R)=2 and op(3,R)=1 and op(4,R)<>1 then RETURN(true): fi: if nops(R)>=5 and op(1,R)=2 and op(2,R)=2 and op(3,R)=1 and op(4,R)=1 and op(5,R)=1 then RETURN(true): fi: if nops(true)>=4 and op(1,R)=2 and op(2,R)=2 and op(3,R)=3 and op(4,R)<>op(3,R) and not(nops(R)>=6 and op(4,R)=op(5,R) and op(5,R)=op(6,R)) then RETURN(true): fi: if nops(R)=2 and op(1,R)=2 and op(2,R)=2 then RETURN(true): fi: if nops(R)=3 and op(1,R)=2 and op(2,R)=2 and op(3,R)>=4 then RETURN(true): fi: if nops(R)>=3 and op(1,R)=2 and op(2,R)=2 and op(3,R)>=4 and op(4,R)<>op(3,R) then RETURN(true): fi: fi: false: end: