back to table · edit · history · where entries came from · files · download
10724 bytes, as of the version from 2026-10-16 09:47 (current). Recorded here, not run.
"""Independent checks for numberdb.org/T460.
sage -python checks.py # compute locally and check
sage -python checks.py --stored # check the database values
Download generate.py and check.py beside this file. --stored needs the owner's
NUMBERDB_API_KEY while T460 is a draft. The campaign wrapper may use
NUMBERDB_KEY_FROM_STDIN=1; that also selects stored checks.
The independent engine implements Mirzakhani's geometric recursion (Do,
arXiv:1103.4674, Theorem 18) in dictionaries of powers of squared lengths.
It shares no intersection arithmetic or partition code with generate.py.
The kernel moments use independently generated Bernoulli numbers, and controls
check F1, F3, V04 and V12 before any general comparison. Published polynomial
coefficients for all eight rows are transcribed from Do, Appendix B.
"""
import numberdb.sage as numberdb
from sage.rings.rational_field import QQ
from sage.rings.polynomial.polynomial_ring_constructor import PolynomialRing
from sage.rings.real_mpfi import RealIntervalField
from fractions import Fraction as F
from functools import lru_cache
from itertools import permutations
from math import comb, factorial
import json
import os
import sys
from pathlib import Path
from check import exactness, agrees_with, names_its_rings, measure
def add(out, key, value):
out[key] = out.get(key,F(0))+value
@lru_cache(None)
def bernoulli(n):
if n == 0:
return F(1)
return -sum((comb(n+1,j)*bernoulli(j) for j in range(n)),F(0))/(n+1)
@lru_cache(None)
def moment(k):
# Coefficient indexed by i means t^(2k-2i)*pi^(2i).
out = {0:F(1,2*k)}
for i in range(1,k+1):
zeta_coefficient = F((-1)**(i+1)*2**(2*i-1),factorial(2*i))*bernoulli(2*i)
out[i] = F(factorial(2*k-1)*(2**(2*i+1)-4),factorial(2*k-2*i))*zeta_coefficient
return out
def double_integral(out, a, b, tail, power, coefficient):
k = a+b+2
beta = F(factorial(2*a+1)*factorial(2*b+1),factorial(2*a+2*b+3))
for i,c in moment(k).items():
add(out,(k-i,)+tail+(power+i,),coefficient*beta*c)
@lru_cache(None)
def geometric(g,n):
# Keys give powers of x_i=L_i^2, followed by P=pi^2.
if g < 0 or 2*g-2+n <= 0:
return {}
if (g,n)==(0,3):
return {(0,0,0,0):F(1)}
if (g,n)==(1,1):
return {(1,0):F(1,48),(0,1):F(1,12)}
if n==0:
out={}
for (a,p),c in geometric(g,1).items():
if a:
add(out,(p+a-1,),c*a*F((-4)**(a-1),g-1))
return out
out={}
for e,c in geometric(g-1,n+1).items():
double_integral(out,e[0],e[1],e[2:-1],e[-1],c)
for mask in range(1 << (n-1)):
I=tuple(i for i in range(n-1) if mask & (1<<i))
J=tuple(i for i in range(n-1) if not mask & (1<<i))
for h in range(g+1):
if 2*h-2+len(I)+1 <= 0 or 2*(g-h)-2+len(J)+1 <= 0:
continue
left=geometric(h,len(I)+1)
right=geometric(g-h,len(J)+1)
for e,c in left.items():
for f,d in right.items():
tail=[0]*(n-1)
for pos,power in zip(I,e[1:-1]): tail[pos]=power
for pos,power in zip(J,f[1:-1]): tail[pos]=power
double_integral(out,e[0],f[0],tuple(tail),e[-1]+f[-1],c*d)
smaller=geometric(g,n-1) if n>1 else {}
for j in range(1,n):
others=tuple(i for i in range(1,n) if i!=j)
for e,c in smaller.items():
for i,f in moment(e[0]+1).items():
r=e[0]+1-i
for b in range(r+1):
powers=[0]*n
powers[0],powers[j]=r-b,b
for pos,power in zip(others,e[1:-1]): powers[pos]=power
add(out,tuple(powers)+(e[-1]+i,),2*comb(2*r,2*b)*c*f)
return {e:c/F(2*(2*e[0]+1)) for e,c in out.items() if c}
def as_polynomial(data,n):
R=PolynomialRing(QQ,names=tuple('L%d'%i for i in range(1,n+1))+('p',))
return R({(tuple(2*x for x in e) if n else 2*e[0]):
QQ(c.numerator)/QQ(c.denominator) for e,c in data.items()})
# (partition of powers of squared lengths, power of pi^2, coefficient).
# These are source data, not evaluations of either recursive implementation.
PUBLISHED={
(0,3):[((),0,'1')],
(0,4):[((1,),0,'1/2'),((),1,'2')],
(0,5):[((2,),0,'1/8'),((1,1),0,'1/2'),((1,),1,'3'),((),2,'10')],
(1,1):[((1,),0,'1/48'),((),1,'1/12')],
(1,2):[((2,),0,'1/192'),((1,1),0,'1/96'),((1,),1,'1/12'),((),2,'1/4')],
(1,3):[((3,),0,'1/1152'),((2,1),0,'1/192'),((1,1,1),0,'1/96'),
((2,),1,'1/24'),((1,1),1,'1/8'),((1,),2,'13/24'),((),3,'14/9')],
(2,0):[((),3,'43/2160')],
(2,1):[((4,),0,'1/442368'),((3,),1,'29/138240'),((2,),2,'139/23040'),
((1,),3,'169/2880'),((),4,'29/192')]
}
def printed_polynomial(g,n):
data={}
for partition,p,c in PUBLISHED[g,n]:
for e in set(permutations(partition+(0,)*(n-len(partition)))):
data[e+(p,)]=F(c)
return as_polynomial(data,n)
def square_coefficients(poly,n):
result={}
for e,c in poly.dict().items():
e=tuple(e) if n else (int(e),)
assert all(int(a)%2==0 for a in e)
result[tuple(int(a)//2 for a in e)]=F(str(c))
return result
def properties(values):
counts={'symmetry_swaps':0,'homogeneous_terms':0,'string_relations':0,
'dilaton_relations':0,'genus_zero_leading_coefficients':0}
data={key:square_coefficients(v,key[1]) for key,v in values.items()}
for (g,n),v in data.items():
dimension=3*g-3+n
for e,c in v.items():
assert c>0 and sum(e)==dimension
counts['homogeneous_terms']+=1
if g==0 and e[-1]==0:
expected=F(factorial(n-3),2**(n-3))
for a in e[:-1]: expected/=factorial(a)**2
assert c==expected
counts['genus_zero_leading_coefficients']+=1
for j in range(n-1):
swapped={e[:j]+(e[j+1],e[j])+e[j+2:]:c for e,c in v.items()}
assert swapped==v
counts['symmetry_swaps']+=1
if (g,n+1) in data:
evaluated,derivative={},{ }
for e,c in data[g,n+1].items():
a,p=e[-2:]
add(evaluated,e[:-2]+(p+a,),c*(-4)**a)
if a: add(derivative,e[:-2]+(p+a-1,),2*a*c*(-4)**(a-1))
rhs={}
for e,c in v.items():
for j in range(n):
f=list(e); f[j]+=1
add(rhs,tuple(f),c/F(2*(e[j]+1)))
strip=lambda x:{e:c for e,c in x.items() if c}
assert strip(evaluated)==strip(rhs),('string',g,n)
assert strip(derivative)=={e:(2*g-2+n)*c for e,c in v.items()},('dilaton',g,n)
counts['string_relations']+=1;counts['dilaton_relations']+=1
return counts
def run_checks(values):
# Test controls first. Failure here invalidates all ensuing comparisons.
assert moment(1)=={0:F(1,2),1:F(2,3)}
assert moment(2)=={0:F(1,4),1:F(2),2:F(28,15)}
assert as_polynomial(geometric(0,4),4)==printed_polynomial(0,4)
assert as_polynomial(geometric(1,2),2)==printed_polynomial(1,2)
expected={(g,n) for g in range(3) for n in range(6) if 0<2*g-2+n<=3}
assert set(values)==expected==set(PUBLISHED)
assert not exactness(values)
assert not agrees_with(values,lambda key:as_polynomial(geometric(*key),key[1]))
assert not agrees_with(values,lambda key:printed_polynomial(*key))
report={'entries':len(values),'independent_geometric_matches':len(values),
'published_polynomial_matches':len(values),'controls':4,
'exactness_complaints':[],'range_complete':True}
report.update(properties(values))
report['measurement']=measure(values)
report['measurement']['longest_at']=str(report['measurement']['longest_at'])
return report
def fetch_values():
tree=numberdb.table('T460')
values={}
for g,rows in tree['Numbers'].items():
for n,entry in rows.items():
text=entry['number'] if isinstance(entry,dict) else entry
ring=PolynomialRing(QQ,names=tuple('L%d'%i for i in range(1,int(n)+1))+('p',))
values[int(g),int(n)]=ring(text)
return values
def corpus_checks(values):
intersections=numberdb.table('T458')['Numbers']
volumes=numberdb.table('T459')['Numbers']
counts={'T458_top_degree_checks':0,'T459_constant_term_checks':0}
RIF=RealIntervalField(512)
assert RIF(1)/RIF(3) > 0 and RIF(1)/RIF(3) < 1
# Decimal spellings denote +/- one final-place unit, not point values.
for (g,n),v in values.items():
coeffs=square_coefficients(v,n)
for e,c in coeffs.items():
if e[-1]==0 and n:
profile=','.join(str(a) for a in sorted(e[:-1]))
stored=intersections[str(g)][str(n)][profile]
stored=stored['number'] if isinstance(stored,dict) else stored
expected=F(stored)
for a in e[:-1]: expected/=2**a*factorial(a)
assert c==expected
counts['T458_top_degree_checks']+=1
dimension=3*g-3+n
c=coeffs[(0,)*n+(dimension,)]
x=(RIF(c.numerator)/RIF(c.denominator))*RIF.pi()**(2*dimension)
text=volumes[str(g)][str(n)]
text=text['number'] if isinstance(text,dict) else text
center=F(text)
if '.' in text:
places=len(text.split('.')[1])
ulp=F(1,10**places)
lo,hi=center-ulp,center+ulp
assert x.lower() >= (RIF(lo.numerator)/RIF(lo.denominator)).upper()
assert x.upper() <= (RIF(hi.numerator)/RIF(hi.denominator)).lower()
assert x.absolute_diameter() < (RIF(ulp.numerator)/RIF(ulp.denominator)/100).lower()
else:
assert n==3 and g==0 and center==1 and x==1
assert not x.is_NaN()
counts['T459_constant_term_checks']+=1
return counts
if __name__=='__main__':
stored='--stored' in sys.argv or os.environ.get('NUMBERDB_KEY_FROM_STDIN')=='1'
if os.environ.get('NUMBERDB_KEY_FROM_STDIN')=='1':
token=sys.stdin.read().strip()
if '=' in token and token.split('=',1)[0].isupper():
token=token.split('=',1)[1].strip().strip('\'"')
numberdb.configure(api_key=token)
if stored:
values=fetch_values()
else:
from generate import WeilPeterssonPolynomials
generator=WeilPeterssonPolynomials()
values={(int(p['g']),int(p['n'])):generator.value(p) for p in generator.enumerate()}
report=run_checks(values)
if stored:
report.update(corpus_checks(values))
report['mode']='stored' if stored else 'local'
print(json.dumps(report,indent=2))