#!/usr/local/bin/maple #-*- maplev -*- # Nathaniel Shar # HW 20 # Experimental Mathematics # It is okay to link to this assignment on the course webpage. Help := proc(): print(`AdjustBFP(a,B,FP,BFP), GFun(A,B,FP,x), GFunT(A,B,FP,T), AllWords(A, n), UA(k), AUA(k), IsImpliedPattern(pat1, pat2, B), MinimalFP(B, FP), Avoids(W, s), SAW(n), SAP(n), GFsawd(d, T)`): end: with(combinat): ################ # From c20.txt # ################ #AdjustBFP(a,B,FP,BFP): inputs a letter a, a symbol B for #blank, a set of forbidden patterns FP, and another set #of forbidden patterns at the beginning #outputs the new BFP for the chopped word assuming #it starts with a AdjustBFP:=proc(a,B,FP,BFP) local NBFP,fp: NBFP:={}: if member([],BFP) then RETURN(FAIL): fi: for fp in FP union BFP do if fp[1]=a or fp[1]=B then NBFP:=NBFP union {fp[2..nops(fp)]}: fi: od: NBFP: end: #GFun(A,B,FP,x): inputs an alphabet A ( a set of symbols or numbers) #and a letter B (NOT in A) and a set of forbidden patterns #FP, and another symbol x (not in A union {B}) and outputs #a rational function in x[A[1]], ..., x[A[nops(A)]] such #the if you mtaylor it the coeff. of #x[A1]^a1*x[A2]^a2*... is the EXACT number of words in #the alphabet A, avoiding the patterns in FP with #a1 occurrences of A1, #a2 occurrences of A2, .... #For example #GFun({1},B,{},x); should output: #1/(1-x[1]); #GFun({1},B,{[1,1]},x); should output: #1+X[1] #weight([1,2,1,1])=x[1]*x[2]*x[1]*x[1] GFun:=proc(A,B,FP,x) local eq,var, F,S, ToDo,BFP,eq1,ABFP,a: #F[BFP]=generating function of all words in the alphabet A #that avoid FP and in addition avoid BFP #S, the set of BFPs that show up #f=1+x[1]*f+x[2]*f+...+x[n]*f #(1-x[1]-...-x[n])*f=1 #f=1/(1-x[1]-...-x[n]) #f=1/(1-n*t) (if x[1]=...=x[n]=t) eq:={}: S:={{}}: var:={}: ToDo:={{}}: while ToDo<>{} do BFP:=ToDo[1]: S:=S union {BFP}: eq1:=F[BFP]-1: var:=var union {F[BFP]}: for a in A do ABFP:=AdjustBFP(a,B,FP,BFP): if not member([],ABFP) then eq1:=eq1-x[a]*F[ABFP]: ToDo:=ToDo union {ABFP}: fi: od: ToDo:=ToDo minus S: eq:=eq union {eq1}: od: var:=solve(eq,var): subs(var,F[{}]): end: #GFunT(A,B,FP,T): inputs an alphabet A ( a set of symbols or numbers) #and a letter B (NOT in A) and a set of forbidden patterns #FP, and another symbol x (not in A union {B}) and outputs #a rational function in x[A[1]], ..., x[A[nops(A)]] such #the if you mtaylor it the coeff. of #T^(length) #For example #GFunT({1},B,{},T); should output: #1/(1-T); #GFun({1},B,{[1,1]},T); should output: #1+T GFunT:=proc(A,B,FP,T) local eq,var, F,S, ToDo,BFP,eq1,ABFP,a: #F[BFP]=generating function of all words in the alphabet A #that avoid FP and in addition avoid BFP #S, the set of BFPs that show up #f=1+x[1]*f+x[2]*f+...+x[n]*f #(1-x[1]-...-x[n])*f=1 #f=1/(1-x[1]-...-x[n]) #f=1/(1-n*t) (if x[1]=...=x[n]=t) eq:={}: S:={{}}: var:={}: ToDo:={{}}: while ToDo<>{} do BFP:=ToDo[1]: S:=S union {BFP}: eq1:=F[BFP]-1: var:=var union {F[BFP]}: for a in A do ABFP:=AdjustBFP(a,B,FP,BFP): if not member([],ABFP) then eq1:=eq1-T*F[ABFP]: ToDo:=ToDo union {ABFP}: fi: od: ToDo:=ToDo minus S: eq:=eq union {eq1}: od: var:=solve(eq,var): subs(var,F[{}]): end: ################################################################################ ############# # Problem 1 # ############# # From HW 18: AllWords := proc(A, n) local w, i: option remember: if n = 1 then: return {seq([i], i in A)}: else: return {seq(seq([op(w), i], i in A), w in AllWords(A, n-1))}: fi: end: UA := proc(k) local R, S, i: R := {}: # i := 1: for S in powerset(AllWords({0,1}, k)) do: if denom(normal(GFunT({0,1}, B, S, x))) = 1 then: R := R union {S}: fi: # print(i): # progress report # i := i + 1: od: return R: end: # UA(2) = {{[0, 0], [0, 1], [1, 1]}, {[0, 0], [1, 0], [1, 1]}, {[0, 0], [0, 1], [1, 0], [1, 1]}} # UA(3): See UA3.txt # UA(4): See UA4.txt ############# # Problem 2 # ############# AUA := proc(k) local R, S, rts, rt, i, good: R := {}: i := 1: for S in powerset(AllWords({0,1}, k)) do: rts := solve(denom(normal(GFunT({0,1}, B, S, x))), x): good := 1: print(S, rts): for rt in {rts} do: if abs(evalf(rt)) < 1 then: good := 0: break: fi: od: if good = 1 then: R := R union {S}: fi: # print(i): i := i + 1: od: return R: end: # AUA(2) = everything except {}, {[0,0]}, {[1,1]} # AUA(3): See AUA3.txt ############# # Problem 3 # ############# IsImpliedPattern := proc(pat1, pat2, B) local beginswith, i: option remember: if pat2 = [] then: return true: fi: if nops(pat1) < nops(pat2) then: return false: fi: beginswith := 1: for i from 1 to nops(pat2) do: if pat1[i] <> pat2[i] and pat2[i] <> B then: beginswith := 0: break: fi: od: if beginswith = 1 then: return true: fi: return IsImpliedPattern(subsop(1=NULL, pat1), pat2, B): end: ############# # Problem 4 # ############# MinimalFP := proc(B, FP) local R, i, j: R := FP: for i in FP do: for j in FP minus {i} do: if IsImpliedPattern(i, j, B) then: R := R minus {i}: break: fi: od: od: return R: end: ############# # Problem 5 # ############# # Assuming W is a self-avoiding walk, checks if W followed by s is # also self-avoiding. Avoids := proc(W, s) local i, j, L, M: L := [op(W), s]: for i from 1 to floor(nops(L)/2) do: M := L[nops(L)-2*i+1..nops(L)]: if add(j, j in select(x->abs(x)=1, M)) = 0 and add(j, j in select(x->abs(x)=2, M)) = 0 then: return false: fi: od: return true: end: SAW := proc(n) local R, W, s: option remember: R := {}: if n = 0 then: return {[]}: fi: for W in SAW(n-1) do: for s in {1,-1,2,-2} do: if Avoids(W, s) then: R := R union {[op(W), s]}: fi: od: od: return R: end: ############# # Problem 6 # ############# # The sequence 1, 4, 12, 36, 100, 284, 780, 2172, 5916, 16268, 44100 # is A001411 in OEIS. ############# # Problem 7 # ############# SAP := proc(n) local L, w, R, s: if n = 0 then return {[]}: fi: if n mod 2 = 1 then return {}: fi: R := {}: for s in {-1, 1, -2, 2} do: L := select(x->nops(select(y->y=s, x)) - nops(select(y->y=-s, x)) = 1 and nops(select(y->y=3-abs(s), x)) = nops(select(y->y=-3+abs(s),x)), SAW(n-1)): for w in L do: if Avoids(w[2..n-1], -s) then: R := R union {[op(w), -s]}: fi: od: od: return R: end: ############# # Problem 8 # ############# GFsawd := proc(d, T) local FP, i: FP := {}: for i from 1 to d do: FP := FP union SAP(2*i): od: return GFunT({1,-1,2,-2}, B, FP, T): end: ############# # Problem 9 # ############# # GFsawd(1, T) = -(T+1)/(3*T-1) # GFsawd(2, T) = -(2*T^2+1+3*T^3+2*T)/(T^3+2*T^2+2*T-1) # GFsawd(3, T) = # -(8*T^8+3*T^7-T^6+8*T^5+3*T^4+4*T^3+2*T^2+2*T+1)/(T+1)/(T^6+T^3-T^2+3*T-1) # GFsawd(4, T) = # -(1+24*T^43+28*T^9+5*T^10+22*T^11-21*T^12-T^13-36*T^14-132*T^15-117*T^16+141*T^35-153*T^17-61*T^18-4*T^19+6*T^3+3*T^2+14*T^5+41*T^7+9*T^6+11*T^42+74*T^21+84 *T^20+32*T^44+24*T^45+8*T^46+49*T^34-80*T^40+6*T^4-117*T^32+47*T^28+2*T^26-117*T^30-94*T^31-82*T^24+24*T^25+61*T^29+16*T^27+101*T^36+52*T^37-16*T^38-72*T^39 +141*T^22-18*T^8-33*T^41+50*T^23-54*T^33+2*T)/(-1-4*T^9+3*T^10-2*T^11+T^12-7*T^13-8*T^14-8*T^15-7*T^16+3*T^35+9*T^17+T^18+12*T^19+2*T^3+T^2+2*T^5-5*T^7-T^6+ T^42+6*T^21+4*T^20+3*T^34+2*T^4+5*T^32-3*T^28+2*T^26-7*T^30-2*T^31-2*T^24+4*T^25-5*T^29+4*T^27-T^36-4*T^37-13*T^22+2*T^8+T^41+2*T^23+2*T^33+2*T) # Only the first is in OEIS.