back to table · edit · history · where entries came from · files · download
7739 bytes, as of the version from 2026-10-29 10:07 (current). Recorded here, not run.
"""Independent checks for upper Khintchine constants, numberdb.org/T606.
sage -python checks.py # local values, certified and numeric checks
sage -python checks.py --stored # repeat on database readback; key for draft
Gaussian quadrature uses mpmath at 160 decimal digits, never Gamma or the
production expression. Known moments 0,2,4 are controls, run first. Its agreement
is corroboration, not certification. Certification comes from Arb and an exact
check that each full ball lies in its written decimal interval. Checks also
include even pairing counts, odd factorials, moment recurrence, monotonicity,
and exact averaging over finite sign vectors. JSON evidence is written to
stdout; progress to stderr. No website write is performed.
"""
import itertools
import json
import os
import sys
from decimal import Decimal
from fractions import Fraction
from pathlib import Path
import numberdb.sage as numberdb
from sage.rings.integer_ring import ZZ
from sage.rings.rational_field import QQ
from sage.rings.real_arb import RealBallField
import sage.version
import mpmath as mp
from numberdb._write import to_text
from numberdb import __version__ as numberdb_version
from generate import UpperKhintchine
from check import exactness, measure, names_its_rings
mp.mp.dps = 160
RB = RealBallField(600)
def rational(text):
q = Fraction(str(text))
return QQ(q.numerator) / QQ(q.denominator)
def decimal_bounds(text):
d = Decimal(text)
midpoint = rational(text)
if '.' not in text and 'e' not in text.lower():
return midpoint, midpoint
unit = QQ(10) ** int(d.as_tuple().exponent)
return midpoint - unit, midpoint + unit
def to_ball(text):
lo, hi = decimal_bounds(text)
return RB(lo).union(RB(hi))
def mp_rational(q):
return mp.mpf(str(q.numerator())) / mp.mpf(str(q.denominator()))
def gaussian_moment(p):
p = mp_rational(p)
return mp.sqrt(2 / mp.pi) * mp.quad(
lambda x: x**p * mp.exp(-x*x/2), [0, 1, 4, 10, 25, mp.inf])
def factorial(m):
ans = ZZ(1)
for k in range(1, m+1):
ans *= k
return ans
def pairing_count(m):
# Choose the partner of the first point, recursively count the rest.
ans = ZZ(1)
for k in range(1, m+1):
ans *= 2*k-1
return ans
def assert_narrow(ball, cutoff=110):
assert ball.is_finite()
assert RB(ball.rad()) < RB(10)**(-cutoff) * abs(ball), str(ball)
def sign_average(vector, p):
total = RB(0)
for signs in itertools.product((-1, 1), repeat=len(vector)):
value = abs(sum(a*s for a,s in zip(vector,signs)))
total += RB(value)**RB(p) if value else RB(0)
return total / RB(2**len(vector))
def main():
if os.environ.get('NUMBERDB_KEY_FROM_STDIN') == '1':
numberdb.configure(api_key=sys.stdin.read().strip())
stored_mode = '--stored' in sys.argv or os.environ.get('CHECK_STORED') == '1'
g = UpperKhintchine()
values = {params['p']: g.value(params, 100) for params in g.enumerate()}
assert not exactness(values)
assert not names_its_rings(Path(__file__).with_name('generate.py'))
texts = {p: to_text(x,digits=100) for p,x in values.items()}
if stored_mode:
doc = numberdb.table('T606')
raw = doc['Numbers']
assert isinstance(raw,dict)
saved = {str(p): str(v['number'] if isinstance(v,dict) else v) for p,v in raw.items()}
assert list(saved) == list(texts), 'Stored enumeration or ordering differs'
texts = saved
controls = {}
for p, known in [(0,1),(2,1),(4,3)]:
value = gaussian_moment(QQ(p))
err = abs(value-known)
assert err < mp.mpf('1e-145')
controls[str(p)] = mp.nstr(err,8)
assert sign_average((1,),QQ(4)) == 1
assert sign_average((1,1),QQ(4)) == 8
assert sign_average((1,-2,3),QQ(2)) == 14
print('Known Gaussian and finite-sign controls passed.',file=sys.stderr,flush=True)
rows = []
max_rel_radius = RB(0)
max_quad_error = mp.mpf(0)
balls = {}
vectors = [(1,), (1,1), (1,2), (1,-2,3), (1,1,1,1)]
for index,(identity,value) in enumerate(values.items(),1):
p = QQ(identity)
text = texts[identity]
lo,hi = decimal_bounds(text)
if p <= 2:
assert text == '1' and value == 1
rel_radius = RB(0)
else:
assert_narrow(value,125)
assert value.lower().exact_rational() >= lo
assert value.upper().exact_rational() <= hi
rel_radius = RB(value.rad()) / abs(RB(value))
ball = to_ball(text)
assert_narrow(ball,98)
balls[p] = ball
moment = gaussian_moment(p)
independent = max(mp.mpf(1), moment**(1/mp_rational(p)))
error = abs(mp.mpf(text)-independent) / independent
max_quad_error = max(max_quad_error,error)
assert mp_rational(lo) <= independent <= mp_rational(hi)
if p.denominator()==1 and p>=2:
n = int(p)
if n%2==0:
expected = RB(pairing_count(n//2))
else:
expected = RB(2)**n * RB(factorial((n-1)//2))**2 / RB.pi()
observed = ball**n if n%2==0 else ball**(2*n)
assert_narrow(expected,90)
assert_narrow(observed,90)
assert observed.overlaps(expected), identity
# Check the inequality directly on five finite coefficient vectors.
for vector in vectors:
norm2 = sum(a*a for a in vector)
if p==2:
assert sign_average(vector,p) == norm2
elif p==1:
# E|S| <= sqrt(E S^2), checked without square-root rounding.
average = sum(abs(sum(a*s for a,s in zip(vector,signs)))
for signs in itertools.product((-1,1),repeat=len(vector)))
exact_average = QQ(average) / QQ(2**len(vector))
assert exact_average**2 <= norm2
else:
ratio = (sign_average(vector,p).log()/RB(p)).exp()/RB(norm2).sqrt()
assert_narrow(ratio,100)
assert ratio < ball
max_rel_radius = max(max_rel_radius,rel_radius)
rows.append({'p':identity,'number':text,'arb_relative_radius':str(rel_radius),
'quadrature_relative_error':mp.nstr(error,10)})
if index%10==0:
print('Independent checks: %d/60'%index,file=sys.stderr,flush=True)
recurrence_count = 0
for p,ball in balls.items():
if p>=2 and p+2 in balls:
left = balls[p+2]**RB(p+2)
right = RB(p+1)*ball**RB(p)
assert_narrow(left,90);assert_narrow(right,90)
assert left.overlaps(right), str(p)
recurrence_count += 1
ordered = sorted(balls)
for a,b in zip(ordered,ordered[1:]):
assert balls[a] < balls[b] if b>2 else balls[a] == balls[b]
print(json.dumps({'table':'T606','mode':'stored' if stored_mode else 'local',
'skill_manifest':'agents/runs/20261029T093155Z-build-skills/manifest.json',
'skill_sha256':'28ecbeeccb0c5e5d6eb541b5a5b0f88aec47b154374d909085116d34f5629199',
'versions':{'sage':sage.version.version,'numberdb':numberdb_version,'mpmath':mp.__version__},
'working_bits':429,'quadrature_digits':160,'entries':len(rows),
'controls_absolute_errors':controls,'finite_sign_checks':len(vectors)*len(rows),
'integer_identity_checks':31,'recurrence_checks':recurrence_count,
'max_arb_relative_radius':str(max_rel_radius),
'max_quadrature_relative_error':mp.nstr(max_quad_error,12),
'measurement_balls':measure(values),'measurement_written':measure(texts),
'rows':rows},indent=2))
if __name__=='__main__':
main()