From: Doron Zeilberger zeilberg AT ekhad.rutgers.edu Date: Mon, 19 Aug 2013 14:23:31 -0400 Subject: source code (Hardin) ez:=proc(): print(`Alp(k), Hkx(k,x), Gkx(k,x)`): print(` IsFollower(w,v), TM(k), Sipur(k,x,L) `): end: #Alp(k): The set of words in {0,1} of length k #with no 00 Alp:=proc(k) local gu,mu,mu1: option remember: if k=0 then RETURN({[]}): fi: if k=1 then RETURN({[0],[1]}): fi: mu:=Alp(k-1): gu:={}: for mu1 in mu do if mu1[nops(mu1)]=0 then gu:=gu union {[op(mu1),1]}: else gu:=gu union {[op(mu1),0],[op(mu1),1]}: fi: od: gu: end: #IsFollower(w,v): Is v a legal follower of w? IsFollower:=proc(w,v) local i,k: k:=nops(w): if nops(v)<>k then RETURN(FAIL): fi: for i from 1 to k do if w[i]=0 and v[i]=0 then RETURN(false): fi: od: for i from 1 to k-1 do if w[i]=0 and v[i+1]=0 then RETURN(false): fi: od: true: end: #TM(k): The transfer matrix of k by n Hardin matrices TM:=proc(k) local i,j,gu,A,A1,w: gu:=convert(Alp(k),list): A:=[]: for i from 1 to nops(gu) do w:=gu[i]: A1:=[]: for j from 1 to nops(gu) do if IsFollower(w,gu[j]) then A1:=[op(A1),1]: else A1:=[op(A1),0]: fi: od: A:=[op(A),A1]: od: A: end: #Gkx(k,x): the generating function, in x, #whose coefficient of x^n is the number of #all k by n Hardin matrices (without restiction #that the top-left entry is 1). See A226444. Gkx:=proc(k,x) local gu,z,eq, eq1,var,i,j: option remember: gu:=Alp(k): var:={seq(z[gu[i]],i=1..nops(gu))}: eq:={}: for i from 1 to nops(gu) do eq1:=z[gu[i]]-x: for j from 1 to nops(gu) do if IsFollower(gu[i],gu[j]) then eq1:=eq1-x*z[gu[j]]: fi: od: eq:=eq union {eq1}: od: var:=solve(eq,var): normal(add(subs(var,z[gu[i]]),i=1..nops(gu))): end: #Hkx(k,x): the generating function, in x, #whose coefficient of x^n is the number of #all k by n Hardin matrices WITH the restiction #that the top-left entry is 1. See A228285. Hkx:=proc(k,x) local gu,z,eq, eq1,var,i,j,mu: option remember: gu:=Alp(k): var:={seq(z[gu[i]],i=1..nops(gu))}: eq:={}: for i from 1 to nops(gu) do eq1:=z[gu[i]]-x: for j from 1 to nops(gu) do if IsFollower(gu[i],gu[j]) then eq1:=eq1-x*z[gu[j]]: fi: od: eq:=eq union {eq1}: od: var:=solve(eq,var): mu:=0: for i from 1 to nops(gu) do if gu[i][1]=0 then mu:=mu+subs(var,z[gu[i]]): fi: od: normal(mu): end: #Sipur(K,x,L): inputs positive integers K and x #and tells the story of the rigorous enumeration #of k by n Hardin matrices, for k from 1 to K #it also prints out #the first L terms of each sequence #and the enumerating sequnce #for square matrices for sizes from 1 to K, #Try: Sipur(6,x); Sipur:=proc(K,x,L) local f, gu,mu,k,i,t0: t0:=time(): print(`Rigorously derived generating functions for`): print(` enumerating k by n Hardin matrices for k from 1`): print(`to `, K): print(``): print(`By Shalosh B. Ekhad `): print(``): print(`Definition: A Hardin matrix has entries from {0,1}`): print(`with no adjacent 1's horizontally, vertically, or`): print(`in the NE-SW direction, and whose top-left entry is 1`): print(``): print(`Definition:For a fixed k, Let f(k)(x) be the `): print(`generating function`): print(`whose coefficient of x^n is the number of n by k`): print(`We have the following (rigorously-derived) facts`): print(`followed by the first`, L, `terms of the enumerating`): print(`sequence (for the sake of Sloane) `): mu:=[]: for k from 1 to K do gu:=Hkx(k,x): lprint(f[k](x)=gu): gu:=taylor(gu,x=0,L+1): print([seq(coeff(gu,x,i),i=1..L)]): mu:=[op(mu),coeff(gu,x,k)]: od: print(`The first`, K, `terms of the enumerating sequence`): print(`for square Hardin matrices is`): print(mu): print(`This took `, time()-t0, `seconds, to generate `): end: #Doron, #My good friend Ron Hardin has been enumerating #matrices of many different classes over the past several years. #His method is to collect a LOT of data, and then to #look for patterns in the numbers. #He just posted the following message to the Sequence #Fans Mailing List. #(start) #The orders of the column (or row) linear recurrences of this s#eem to form #the Fibonacci series. # #T(n,k)=Number of nXk binary arrays with top left value 1 and no two ones #adjacent horizontally, vertically or nw-se diagonally #Table starts #..1...1.....2......3........5.........8..........13......... #..21.............34 #..1...1.....3......6.......13........28..........60......... #.129............277 #..2...3....13.....35......112.......337........1034......... #3154...........9637 #..3...6....35....133......587......2448.......10414........ #44024.........186414 #..5..13...112....587.....3631.....21166......126119....... #745178........4416695 #..8..28...337...2448....21166....172082.....1428523..... #11771298.......97268701 #.13..60..1034..10414...126119...1428523....16566199.... #190540884.....2197847780 #.21.129..3154..44024...745178..11771298...190540884... #3057290265....49208639399 #.34.277..9637.186414..4416695..97268701..2197847780.. #49208639399..1105411581741 #.55.595.29431.789100.26150120.802886174.25325358687. #791176762937.24801939723742 #Column 1 is A000045 #Column 2 is A002478(n-1) #Empirical for column k: # k=1: a(n)=a(n-1)+a(n-2) #k=2: a(n)=a(n-1)+2*a(n-2)+a(n-3) #k=3: a(n)=a(n-1)+5*a(n-2)+4*a(n-3)-a(n-5) #k=4: a(n)=a(n-1)+10*a(n-2)+15*a(n-3)+4*a(n-4)-6*a(n-5)-a(n-6)+3* #a(n-7)-a(n-8) #k=5: [order 13] #k=6: [order 21] #k=7: [order 34] #k=8: [order 55] #k=9: [order 89] #(end)