Trees

{{{id=48| class tree: # tree represents a rooted tree def __init__(self, branches): # construct a rooted tree by its branches m1=0 self.m0=1 self.order=1 for branch in branches: m1=m1*branch.m+self.m0*branch.m0 self.m0*=branch.m self.order+=branch.order branch_types=union(map(lambda x:x.tree_type,branches)) if branch_types in ([], ["B"]): self.tree_type="A" self.m=self.m0+m1 elif branch_types==["A"]: self.tree_type="B" self.m=m1 else: raise NotImplementedError("branches have types %s" % branch_types) def __repr__(self): return("Tree of order %s with (m,m0)=(%s,%s)" % (self.order, self.m, self.m0)) def rho(self): return(self.m0/self.m) def abtree(branches): # shorthand notation return(tree([tree(branches)])) L=tree([]) F=abtree([L,L]) # L and F as in the paper def C(length, T): assert length in ZZ assert length>= 0 if length>0: return(abtree([L,F,C(length-1,T)])) return(T) /// }}}

Replacement Rules

Important Rooted Trees (Figure 8)

{{{id=65| A_3=abtree([L]) A_6=abtree([L,L,L,L]) A_7=abtree([L,F]) A_10=abtree([L,A_7]) A_14=abtree([F,F,F]) A_24=abtree([F,F,A_14]) /// }}} {{{id=63| assert (L.order,L.rho())==(1,1) assert (A_3.order,A_3.rho())==(3,1/2) assert (F.order,F.rho())==(4,2/3) assert (A_6.order,A_6.rho())==(6,4/5) assert (A_7.order,A_7.rho())==(7,5/8) assert (A_10.order,A_10.rho())==(10,13/21) assert (A_14.order,A_14.rho())==(14,2/3) assert (A_24.order,A_24.rho())==(24,2/3) /// }}}

Replacements for some alpha (Table 3)

{{{id=61| def check_replacement_alpha(T1,T2, alphamin, alphamax, leftopen): # assert that T1.m+T1.m0*alphaT1.m0 and T1.m+alphamin*T1.m0==T2.m+alphamin*T2.m0 else: raise NotImplementedError("alphamin=%s, alphamax=%s, leftopen=%s" % (alphamin, alphamax,leftopen)) /// }}} {{{id=69| check_replacement_alpha(tree([A_3]), tree([L, L, L]), 0,2,False) check_replacement_alpha(tree([L, A_3]), tree([L, L, L, L]), 0,1,False) check_replacement_alpha(tree([A_3, A_3]), tree([L, abtree([L, L, L])]), 0,+Infinity,False) check_replacement_alpha(tree([A_3, A_3, A_3]), tree([L, F, F]), 0,+Infinity,False) check_replacement_alpha(abtree([F, F, A_24]), abtree([L, C(1,L), abtree([L, C(1,L), abtree([L, L, C(1,L)])])]), 50/2473,+Infinity,True) check_replacement_alpha(abtree([F, A_14, A_14]), abtree([L, C(1,L), abtree([L, C(1,L), abtree([L, L, C(1,L)])])]), 50/2473,+Infinity,True) check_replacement_alpha(abtree([F, A_14, A_24]), C(1,abtree([L, C(1,L), abtree([L, C(1,L), C(2,L)])])), 0,+Infinity,False) check_replacement_alpha(abtree([F, A_24, A_24]), C(2,abtree([L, C(2,L), C(3,L)])), 0,+Infinity,False) check_replacement_alpha(abtree([A_14, A_14, A_14]), C(1,abtree([L, C(1,L), abtree([L, C(1,L), C(2,L)])])), 0,+Infinity,False) check_replacement_alpha(abtree([A_14, A_14, A_24]), C(2,abtree([L, C(2,L), C(3,L)])), 0,+Infinity,False) check_replacement_alpha(abtree([A_14, A_24, A_24]), C(9,L), 0,+Infinity,False) check_replacement_alpha(abtree([A_24, A_24, A_24]), C(10,F), 0,+Infinity,False) check_replacement_alpha(tree([L, A_10]), tree([L, L, abtree([L, L, abtree([L, L, L])])]), 0,7/6,False) check_replacement_alpha(tree([L, A_14]), tree([L, L, abtree([L, L, abtree([L, L, abtree([L, L, L])])])]), 0,18/25,False) check_replacement_alpha(tree([L, A_24]), tree([L, L, abtree([L, L, abtree([L, C(1,L), C(1,L)])])]), 0,295/318,False) check_replacement_alpha(tree([A_10, A_7, A_7]), tree([L, abtree([L, L, abtree([L, C(1,L), C(1,L)])])]), 0,+Infinity,False) check_replacement_alpha(tree([A_10, A_10, A_7]), tree([F, C(1,L), C(2,L)]), 0,+Infinity,False) check_replacement_alpha(tree([A_10, A_10, A_10]), tree([F, C(1,L), C(2,F)]), 0,+Infinity,False) check_replacement_alpha(tree([A_7, F, F]), tree([L, L, abtree([L, L, abtree([L, L, abtree([L, L, L])])])]), 0,3/4,False) check_replacement_alpha(tree([A_7, A_7, F]), tree([L, L, abtree([L, L, abtree([L, L, C(1,L)])])]), 0,50/39,False) check_replacement_alpha(tree([A_7, A_7, A_7]), tree([L, abtree([L, L, abtree([L, L, abtree([L, L, C(1,L)])])])]), 0,+Infinity,False) check_replacement_alpha(tree([L, L, L, abtree([L, L, L])]), tree([L, A_7]), 1/2,+Infinity,True) check_replacement_alpha(tree([L, L, abtree([L, L, abtree([L, L, abtree([L, L, abtree([L, L, L])])])])]), tree([L, F, A_14]), 2/17,+Infinity,True) check_replacement_alpha(tree([L, F, C(1,A_14)]), tree([L, F, abtree([F, F, C(1,F)])]), 1/3,+Infinity,True) check_replacement_alpha(tree([L, F, C(1,A_24)]), tree([L, F, abtree([F, F, abtree([F, F, C(1,F)])])]), 1/3,+Infinity,True) /// }}}

Replacements for alpha=0 (Table 2)

{{{id=62| def check_replacement(T1,T2): assert T1.order==T2.order assert T1.mSetup for chains (cf. Lemma 6.2) {{{id=2| var('k') A=matrix([[8,3],[5,3]]) NF=QQ.extension(A.charpoly(),'alpha') alpha=NF.gen() # alpha here corresponds to lambda in the paper, as lambda is a # reserved word in python alphabar=9/alpha # alphabar is the conjugate of alpha assert A.charpoly()(alpha)==0 assert A.charpoly()(alphabar)==0 assert alpha.complex_embedding(i=1)>10 assert alphabar.complex_embedding(i=1)<1 MNF=MatrixSpace(NF,2,2) D,T=MNF(A).eigenmatrix_right() assert A==T*D*T.inverse() assert D==Matrix([[alpha,0],[0,alphabar]]) NFP.=PolynomialRing(NF,2) # alpha_k=alpha^k # alphabar_k=alphabar^k NFRF=NFP.fraction_field() def val(expression,kk): # substitute k=kk in the expression, in particular substitute alpha_k and alphabar_k return(expression.subs({alpha_k:alpha^kk,alphabar_k: (9/alpha)^kk})) D_k=Matrix([[alpha_k,0],[0,alphabar_k]]) # D_k = D^k assert val(D_k,1)==D def coefficients_k(exponent): # determine pair (c0,c1) of integers such that exponent=c0+c1*k if exponent in ZZ: return((exponent,0)) c1=exponent.coefficient(k,1) c0=exponent.coefficient(k,0) assert exponent==c1*k+c0, "%s is not a affine linear function of k" % exponent assert c1 in ZZ, "Coefficient of k in %s is not an integer" % exponent assert c0 in ZZ, "Constant term in %s is not an integer" % exponent return(c0,c1) def alphapower(exponent): # compute alpha^exponent where exponent is an affine linear function in k (c0,c1)=coefficients_k(exponent) return( alpha_k^c1*alpha^c0 ) def Apower(exponent): # compute A^exponent where exponent is an affine linear function in k (c0,c1)=coefficients_k(exponent) return( T*D_k^c1*T.inverse()*A^c0 ) for a in range(-1,1): for b in range(-1,1): for kk in range(3): assert val(Apower(a*k+b),kk)==A^(a*kk+b) /// }}} {{{id=3| class chain(tree): # chain represents a rooted tree which arises by a chain operation def __init__(self, length, T): # construct C^length T for some rooted tree T self.m,self.m0= Apower(length)*vector([T.m,T.m0]) self.order=T.order+7*length assert T.tree_type=="A" self.tree_type="A" def CL(length): # C^length L return(chain(length,L)) def CF(length): # C^length F return(chain(length,F)) /// }}}

Lemma 6.7

Preparations

{{{id=33| def is_positive(nfelement): # returns nfelement>0 for a nfelement in NF if nfelement in QQ: return(QQ(nfelement)>0) assert nfelement.parent()==NF [c0,c1]=nfelement.list() compare=-c0/c1 if c1>0: return(not(compare>10 and A.characteristic_polynomial()(compare)>0)) else: return(compare>10 and A.characteristic_polynomial()(compare)>0) def is_negative(nfelement): # returns nfelement<0 for a nfelement in NF return(is_positive(-nfelement)) def absNF(nfelement): # returns abs(nfelement) if is_negative(nfelement): return(-nfelement) else: return(nfelement) def minNF(l): # returns min(l) for a list l of elements of NF if len(l)==1: return(l[0]) else: a=l[0] b=minNF(l[1:]) if is_positive(a-b): return(b) else: return(a) def sign_of_k_polynomial(what): # returns (sign,kmin) if sign(what)==sign for k>=kmin; here what in NFP assert what.parent()==NFP assert what.is_homogeneous() degree=what.degree() degree_alpha_k=what.degree(alpha_k) leading_coefficient=what.monomial_coefficient( alpha_k^degree_alpha_k*alphabar_k^(degree-degree_alpha_k)) kk=0 while True: extra=add([absNF( what.monomial_coefficient(alpha_k^i*alphabar_k^(degree-i)))* (alpha^(kk*(i-degree_alpha_k))*alphabar^(kk*(degree_alpha_k-i))) for i in range(degree_alpha_k)]) if is_positive(absNF(leading_coefficient)-extra): break kk+=1 if is_positive(leading_coefficient): sign=1 elif is_negative(leading_coefficient): sign=-1 return(sign,kk) def sign_of_fraction(what): # returns (sign,kmin) if sign(what)==sign for k>=kmin; # here what in Quotient Field of NFP signs=[sign_of_k_polynomial(what.numerator()), sign_of_k_polynomial(what.denominator())] return(signs[0][0]*signs[1][0],max(signs[0][1],signs[1][1])) def k_limit(rational_function): # limit of a rational function in alpha_k, alphabar_k, alpha for k to infinity num=rational_function.numerator() den=rational_function.denominator() assert num.parent()==NFP assert den.parent()==NFP assert num.is_homogeneous() assert den.is_homogeneous() assert num.degree(alpha_k)==den.degree(alpha_k) d=num.degree(alpha_k) return(num.monomial_coefficient(alpha_k^d)/den.monomial_coefficient(alpha_k^d)) def find_minimum(sequence): # For a sequence in NFRF, find its minimum over all k>=0. # The sequence has to be monotonic for almost all k and convergent. if sequence in NF: return(sequence) assert sequence in NFRF lim=k_limit(sequence) sign_of_difference,k_stable_sign=sign_of_fraction( sequence.subs({alpha_k:alpha*alpha_k,alphabar_k:alphabar*alphabar_k})-sequence) candidates=[val(sequence,kk) for kk in range(k_stable_sign+1)] candidates.append(lim) return(minNF(candidates)) rhomax=ceil((sqrt(3)-1)*10000)/10000 assert rhomax>sqrt(3)-1 def replace_subtree(old,new): # Let old and new be rooted trees, possibly depending on k. # Returns a pair (q,new.order-old.order) such that # m(Tnew)/m(Told)>= q for all trees Told fulfilling the LC # containing old and the tree Tnew arises from Told # by replacing one occurrence of old in Told by Tnew rhomin=2/3 min_rhomin=find_minimum((new.m+new.m0*rhomin)/(old.m+old.m0*rhomin)) min_rhomax=find_minimum((new.m+new.m0*rhomax)/(old.m+old.m0*rhomax)) return((minNF([min_rhomin, min_rhomax]),new.order-old.order)) def minNFCount(list): # Given a list of pairs (q_i, d_i) (for q_i in NF), determine (min(q_i), union(d_i)) Qmin=minNF(map(lambda x:x[0],list)) if Qmin in QQ: Qmin_rounded=floor(QQ(Qmin)*10000)/10000 else: Qmin_rounded=floor(Qmin.complex_embedding(i=1)*10000)/10000 assert is_positive(Qmin-Qmin_rounded) return(Qmin, Qmin_rounded.n(),union(map(lambda x:x[1],list))) def combine_replacements(list): return(is_positive(prod(map(lambda x:x[0],list))-1), prod(map(lambda x:x[1],list)), add(map(lambda x:x[2][0],list))) /// }}}

Verification of the items

{{{id=34| r1=minNFCount([replace_subtree( tree([CF(2*k+floor((s)/3)),CF(2*k+floor((s+1)/3)), CF(2*k+floor((s+2)/3))]), tree([CL(3*k+floor((s+2)/2)),CL(3*k+floor((s+3) /2)),L])) for s in range(6)]) r1 /// }}} {{{id=38| r2=minNFCount([ replace_subtree( tree([CF(2*k+floor((s)/3)),CF(2*k+floor((s+1)/3)), CF(2*k+floor((s+2)/3))]), tree([CL(3*k+floor((s+1)/2)),CL(3*k+floor((s+2)/2)),L])) for s in range(6)]) r2 /// }}} {{{id=39| r3=minNFCount([replace_subtree( tree([CF(k+floor(s/3)),CF(k+floor((s+1)/3)),CF(k+floor((s+2)/3))]), tree([CF(3*k+s),F,L])) for s in range(3)]) r3 /// }}} {{{id=40| r4=minNFCount([replace_subtree( tree([CF(k+floor((jj+1)/3)),CF(k+floor((jj+2)/3)), chain(k, abtree([CF(k+floor((ii+1)/4)),CF(k+floor((ii+2)/4)),CF(k+floor((ii+3)/4))]))]), tree([L,F,chain(2*k+ii+jj,abtree([CL(2*k+1),CL(2*k+1),L]))])) for ii in range(4) for jj in range(3)]); r4 /// }}} {{{id=41| r5=minNFCount([replace_subtree( tree([CL(k+1),CL(k+j),L]), tree([L,F,CF(2*k+j)])) for j in [0,1]]) r5 /// }}} {{{id=42| r6=minNFCount([replace_subtree( tree([CL(3*k+floor(s/2)),CL(3*k+floor((s+1)/2)),L]), tree([CF(2*k+floor(s/3)),CF(2*k+floor((s+1)/3)),CF(2*k+floor((s-1)/3))])) for s in range(1,7)]) r6 /// }}} {{{id=43| r6a=minNFCount([replace_subtree( tree([CL(3*k+floor(s/2)),CL(3*k+floor((s+1)/2)),L]), tree([CF(2*k+floor(s/3)),CF(2*k+floor((s+1)/3)),CF(2*k+floor((s-1)/3))])) for s in range(2,8)]) r6a /// }}} {{{id=44| r7=minNFCount([replace_subtree( tree([CL(k+2+t),L,chain(k+1, abtree([CL(k+2+floor((s+1)/3)),CL(k+2+floor((s+2)/3)),L]))]), tree([F,L,chain(k+floor((s+t+2)/4), abtree([CF(k+1+floor((s+t+3)/4)),CF(k+1+floor((s+t+4)/4)), CF(k+1+floor((s+t+5)/4))]))])) for s in range(3) for t in range(2)]) r7 /// }}} {{{id=45| r8=minNFCount([replace_subtree( tree([CL(t),L,abtree([CL(floor(s/2)),CL(floor((s+1)/2)),L])]), tree([CF(floor((s+t-1)/3)),CF(floor((s+t)/3)),CF(floor((s+t+1)/3))])) for s in range(1,5) for t in range (0,3)]) r8 /// }}} {{{id=46| r9=minNFCount([replace_subtree( tree([CL(t),L,abtree([CL(floor(s/2)),CL(floor((s+1)/2)),L])]), tree([L,F,CF(s+t)])) for s in range(1,5) for t in range(3)]); r9 /// }}}

Applications of Lemma 6.7

Lemma 6.10

{{{id=47| combine_replacements([r1,r2,r3]) /// }}}

Lemma 6.11

{{{id=50| combine_replacements([r2,r4]) /// }}}

Lemma 6.14

{{{id=51| combine_replacements([r5,r5,r6]) /// }}}

Lemma 6.16

{{{id=52| combine_replacements([r6a,r7]) /// }}} {{{id=53| combine_replacements([r9,r7]) /// }}} {{{id=54| combine_replacements([r8,r9]) /// }}}

Optimal Trees

{{{id=4| optimal_tree_symbolic=[None for kk in range(7)] # optimal_tree_symbolic[case] represents the "generic" case of T_n^* # for order n>=37 with case = n mod 7 # optimal_tree_symbolic[case] is a list containing various subcases, # as we resolve expressions of the shape floor( (r*k+j)/r ) # explicitly before storing the result. optimal_tree_symbolic[1]=[CL(k)] optimal_tree_symbolic[2]=[tree([CL(0),CL(k+1+floor((j+1)/5)),CL(k+1+floor((j+3)/5)), chain(k+floor(j/5), abtree([CL(0),CL(k+1+floor((j+2)/5)),CL(k+1+floor((j+4)/5))]))]) for j in range(5)] optimal_tree_symbolic[3]=[tree([CF(k+floor((j+0)/4)),CF(k+floor((j+1)/4)), CF(k+floor((j+2)/4)),CF(k+floor((j+3)/4))]) for j in range(4)] optimal_tree_symbolic[4]=[CF(k)] optimal_tree_symbolic[5]=[tree([CL(0),CL(k+floor((j+0)/3)),CL(k+floor((j+1)/3)), CL(k+floor((j+2)/3))]) for j in range(3)] optimal_tree_symbolic[6]=[tree([CF(k+floor((j+1)/7)),CF(k+floor((j+3)/7)), CF(k+floor((j+5)/7)), chain(k+floor(j/7), abtree([CF(k+floor((j+2)/7)),CF(k+floor((j+4)/7)),CF(k+floor((j+6)/7))]))]) for j in range(7)] optimal_tree_symbolic[0]=[tree([CL(0),CL(0),CF(k)])] /// }}} {{{id=12| def optimal_m(n): # determine m(T_n^*) for n>=37 for a specific n case=n % 7 # compute relevant subcase and corresponding value of k offset=optimal_tree_symbolic[case][0].order.subs({k:0}) num_subcases=len(optimal_tree_symbolic[case]) subcase=int(((n-offset)/7)) % num_subcases kk = (n-offset-7*subcase)/(7*num_subcases) assert kk>=0, "T_%d^* is not a generic case " % n assert n==optimal_tree_symbolic[case][subcase].order.subs({k:kk}) return(ZZ(val(optimal_tree_symbolic[case][subcase].m,kk))) /// }}}

Lemma 6.16, case k_0=k_5=0, s=1

{{{id=59| open_cases_6_16=[tree([CL(0),CL(floor(s/2)),CL(floor((s+1)/2)), abtree([CL(0),CL(t), abtree([CL(0),CL(floor(u/2)),CL(floor((u+1)/2))])])]) for t in range(3) for s in range(1,5) for u in range(1,5)] for t in open_cases_6_16: assert optimal_m(t.order)>t.m /// }}}

Asymptotics of m(T_n^*)

{{{id=29| def asymptotic_coefficient(case): # returns c_case such that m(T_n^*) ~ c_case lamba^floor(n/7) for all n # with case = n mod 7. result=union([ k_limit(T.m/alphapower((T.order()-case)/7)) for T in optimal_tree_symbolic[case] ]) assert len(result)==1 return(result[0]) /// }}} {{{id=30| def numeric_asymptotic_coefficient(case): # returns c_case_numeric such that m(T_n^*) ~ c_case lamba^(n/7) # for all n with case = n mod 7. return(asymptotic_coefficient(case).complex_embedding(i=1)/ (alpha.complex_embedding(i=1))^(case/7)) /// }}} {{{id=31| numeric_asymptotic_coefficients=map(numeric_asymptotic_coefficient, range(7)) def optimal_m_approx(n): return(numeric_asymptotic_coefficients[ n % 7]*alpha.complex_embedding(i=1)^(n/7)) /// }}} {{{id=32| approx_ratios=[ optimal_m_approx(nn)/optimal_m(nn) for nn in range(1000,1100)] (min(approx_ratios),max(approx_ratios)) /// }}} {{{id=23| [asymptotic_coefficient(jj) for jj in range(7)] /// }}} {{{id=25| [numeric_asymptotic_coefficient(jj) for jj in range(7)] /// }}}