back to table · edit · history · where entries came from · files · download
7795 bytes, as of the version from 2026-10-29 19:37 (current). Recorded here, not run.
"""Independent checks of Löbell growth rates -- numberdb.org/T612.
sage -python checks.py # all local values and identities
sage -python checks.py --stored # also the written decimal intervals
The optional key for a private draft is NUMBERDB_API_KEY. For the campaign,
NUMBERDB_KEY_FROM_STDIN=1 reads it from stdin; NUMBERDB_CHECK_STORED=1
selects the read-back check. No writes are made.
The dual graph is constructed from the two pentagonal flowers (Inoue §4),
then its cliques give finite special subgroups. Steinberg's formula constructs
the denominator without using the generator's radical formula. Its reciprocal
polynomial is bisected in exact rational arithmetic. Short element counts use
only involutions and commutations, independently of that rational series.
"""
import os
import sys
import json
import itertools
from decimal import Decimal
from fractions import Fraction
import numberdb.sage as numberdb
from numberdb._write import to_text
from sage.rings.integer_ring import ZZ
from sage.rings.rational_field import QQ
from sage.rings.polynomial.polynomial_ring_constructor import PolynomialRing
from sage.rings.real_arb import RealBallField
from generate import LobellGrowth
def graph(n):
# 0 and 1 are the bases. Upper faces: 2+i; lower: 2+n+i.
adj = [set() for _ in range(2*n+2)]
def edge(a, b):
adj[a].add(b); adj[b].add(a)
for i in range(n):
u, v = 2+i, 2+n+i
edge(0, u); edge(1, v)
edge(u, 2+(i+1)%n); edge(v, 2+n+(i+1)%n)
edge(u, v); edge(u, 2+n+(i-1)%n)
return adj
def clique_counts(adj):
counts = [1]
cliques = [()]
while cliques:
extended = []
for clique in cliques:
start = clique[-1]+1 if clique else 0
for v in range(start, len(adj)):
if all(v in adj[u] for u in clique):
extended.append(clique+(v,))
if extended:
counts.append(len(extended))
cliques = extended
return counts
def canonical(word, adj):
# Closure under all adjacent commutations and ss -> empty. This exhaustive
# algorithm is only used on words of length at most 3.
seen = {word}; pending = [word]
while pending:
w = pending.pop()
for i in range(len(w)-1):
if w[i] == w[i+1]:
changed = w[:i]+w[i+2:]
elif w[i+1] in adj[w[i]]:
changed = w[:i]+(w[i+1],w[i])+w[i+2:]
else:
continue
if changed not in seen:
seen.add(changed); pending.append(changed)
return min(seen, key=lambda w:(len(w),w))
def short_counts(adj):
layer = {()}; counts = [1]
for length in range(1,4):
nxt = set()
for w in layer:
for s in range(len(adj)):
reduced = canonical(w+(s,), adj)
if len(reduced) == length:
nxt.add(reduced)
counts.append(len(nxt)); layer = nxt
return counts
def decimal_interval(text):
d = Decimal(text)
c = Fraction(d); e = d.as_tuple().exponent
unit = Fraction(10**e) if e >= 0 else Fraction(1,10**(-e))
return QQ(c.numerator)/c.denominator, QQ(unit.numerator)/unit.denominator
def bisect(poly, lo, hi, steps=440):
lo,hi=QQ(lo),QQ(hi)
assert poly(lo)<0<poly(hi)
for _ in range(steps):
middle=(lo+hi)/2
if poly(middle)<0: lo=middle
else: hi=middle
assert poly(lo)<0<poly(hi)
return lo,hi
def unit_power(n, R):
rad=ZZ(n*(n-2)); d=rad.squarefree_part(); b=ZZ(rad//d).isqrt()
assert d*b*b == rad and d>1
# All units 1<epsilon<=tau have positive trace a <= trace(tau)+1;
# a^2-d*b^2=+/-4, with the maximal-order parity conditions.
candidates=[]
for a in range(1,2*n+1):
for c in range(1,2*n+1):
if d%4 == 1:
if (a-c)%2: continue
elif a%2 or c%2: continue
if abs(a*a-d*c*c) == 4:
val=(R(a)+R(d).sqrt()*c)/2
if val>1: candidates.append((val,QQ(a)/2,QQ(c)/2))
assert candidates
_,u,v=min(candidates,key=lambda row:row[0].mid())
x,y=QQ(1),QQ(0)
for power in range(1,100):
x,y=x*u+d*y*v,x*v+y*u
if (x,y)==(n-1,b):
D=d if d%4==1 else 4*d
return int(D),power,(u,v,int(d))
raise AssertionError('No exact fundamental-unit power found')
def main(stored=False):
# Controls before any use of the counting or bracketing procedures.
assert short_counts([{1,2},{0,2},{0,1}]) == [1,3,3,1]
assert short_counts([set(),set(),set()]) == [1,3,6,12]
P=PolynomialRing(QQ,'t'); t=P.gen()
lo,hi=bisect(t*t-2,1,2,80)
assert lo*lo<2<hi*hi and hi-lo==QQ(1)/2**80
print('Controls passed: (Z/2)^3, three free involutions, sqrt(2) bracket.',flush=True)
R=RealBallField(500); g=LobellGrowth(); rows=[]; values={}
doc=numberdb.table('T612') if stored else None
regulators=numberdb.table('T131')['Numbers'] if stored else None
if stored: assert list(doc['Numbers'])==[str(n) for n in range(5,31)]
for params in g.enumerate():
n=int(params['n']); result=g.value(params,100)
ball=result['number'] if isinstance(result,dict) else result
adj=graph(n); c=clique_counts(adj)
assert c==[1,2*n+2,6*n,4*n]
D=sum((-1)**k*c[k]*t**k*(1+t)**(3-k) for k in range(4))
expected=(1-t)*(1-(2*n-2)*t+t*t)
assert D==expected
# x^3 D(1/x) is assembled by reversing coefficients.
reciprocal=P(list(reversed(D.list())))
lo,hi=bisect(reciprocal,2*n-3,2*n-2)
assert ball.is_finite() and R(ball.rad())/R(ball.mid()) < R(10)**-117
assert R(ball).upper()>=R(lo).lower() and R(ball).lower()<=R(hi).upper()
q=ball*ball-(2*n-2)*ball+1
assert q.is_finite() and q.contains_zero() and R(q.rad())<R(10)**-112
assert (ball+1/ball-(2*n-2)).contains_zero()
assert (ball.log()-R(n-1).arccosh()).contains_zero()
assert ZZ(n*(n-2)).is_square() is False
coeff=[]
for k in range(4):
coeff.append(int(((1+t)**3)[k])-sum(D[j]*coeff[k-j] for j in range(1,min(3,k)+1)))
counts=short_counts(adj)
assert coeff==counts
disc,m,unit=unit_power(n,R)
unit_ball=R(unit[0])+R(unit[1])*R(unit[2]).sqrt()
assert (ball-unit_ball**m).contains_zero()
text=to_text(ball,digits=100)
if stored:
entry=doc['Numbers'][str(n)]; text=entry['number'] if isinstance(entry, dict) else entry
if str(disc) in regulators:
reg=regulators[str(disc)]; rc,ru=decimal_interval(reg['number'] if isinstance(reg, dict) else reg)
log_value=R(ball).log()/m
assert R(rc-ru)<log_value<R(rc+ru)
center,ulp=decimal_interval(text)
left,right=center-ulp,center+ulp
assert reciprocal(left)<0<reciprocal(right)
# Containment, not merely overlap, for the serialized ball enclosure.
assert R(left)<R(ball) < R(right)
assert left<lo<hi<right
values[n]=text
rows.append({'n':n,'cliques':c,'short_counts':counts,'D':disc,'unit_power':m,
'relative_accuracy_bits':max(k for k in range(300, 430) if R(ball.rad())/R(ball.mid()) < R(2)**(-k)),'text':text})
print('n=%d: incidence, short counts, rational bracket, unit power, serialization OK'%n,flush=True)
print(json.dumps({'stored':stored,'rows':rows,'minimum_accuracy_bits':min(r['relative_accuracy_bits'] for r in rows)},indent=2))
if __name__=='__main__':
if os.environ.get('NUMBERDB_KEY_FROM_STDIN')=='1':
numberdb.configure(api_key=sys.stdin.read().strip())
main('--stored' in sys.argv or os.environ.get('NUMBERDB_CHECK_STORED')=='1')