back to table · edit · history · where entries came from · files · download
8991 bytes, as of the version from 2026-11-01 09:03 (current). Recorded here, not run.
"""Independent checks for Báez-Duarte coefficients -- numberdb.org/T644.
Run `sage -python checks.py` before filling, or `sage -python checks.py --stored`
after filling (a private draft requires NUMBERDB_API_KEY). Outputs JSON evidence.
The first method uses exact Bernoulli numbers built from their generating
function, an independently enclosed Machin formula for pi, and polynomial
Horner evaluation. It never calls zeta. Every coefficient is checked.
The stored mode verifies enclosure containment, not rounded-midpoint residuals.
Binomial inversion is ill-conditioned at high orders: its stored-value check
has at least 60 relative bits, not 100 decimal digits. This is supplementary
to the independent enclosure check of every 100-digit value.
"""
import json
import sys
import os
import time
from decimal import Decimal
from fractions import Fraction
from math import comb, factorial
import numberdb.sage as numberdb
from sage.rings.rational_field import QQ
from sage.rings.real_arb import RealBallField
from sage.rings.real_mpfi import RealIntervalField
from sage.env import SAGE_VERSION
from numberdb._write import to_text
from numberdb._compare import _claim, digits_of
from numberdb._generate import _stored_entries
from generate import BaezDuarte, MAX_K, coefficient_balls, configure_key_from_stdin
from check import exactness, measure, names_its_rings
def bernoulli_recurrence(limit):
# (exp(t)-1) sum B_n*t^n/n! = t, including B_1=-1/2.
result = [QQ(1)]
for n in range(1, limit + 1):
if n > 1 and n % 2:
result.append(QQ(0))
else:
result.append(-sum((QQ(comb(n+1,j))*result[j]
for j in range(n)), QQ(0)) / QQ(n+1))
return result
def arctangent_bounds(q, bits):
# Alternating decreasing series; the exact next term bounds the error.
n = 0
total = QQ(0)
while True:
total += QQ((-1)**n) / QQ((2*n+1)*q**(2*n+1))
n += 1
term = QQ((-1)**n) / QQ((2*n+1)*q**(2*n+1))
if abs(term) < QQ(1) / QQ(2**(bits+20)):
return min(total,total+term), max(total,total+term), n
def independent(limit):
precision = numberdb.bits(150, losing=limit+128)
field = RealBallField(precision)
intervals = RealIntervalField(precision)
tangent=QQ(1)/5
for _ in range(2):
tangent=2*tangent/(1-tangent*tangent)
assert (tangent-QQ(1)/239)/(1+tangent/239) == 1
a,b,n5=arctangent_bounds(5,precision)
c,d,n239=arctangent_bounds(239,precision)
pi = field(16*intervals(a,b)-4*intervals(c,d))
assert pi.is_finite() and pi.overlaps(field.pi())
assert pi > 3 and pi < 4
bern = bernoulli_recurrence(2*limit+2)
controls = {0:QQ(1),1:-QQ(1)/2,2:QQ(1)/6,4:-QQ(1)/30,
6:QQ(1)/42,12:-QQ(691)/2730}
assert all(bern[n] == value for n,value in controls.items())
assert not exactness(bern)
# Euler's even-zeta formula gives f_j = a_j * t^(j+1), t = pi^(-2).
a = [QQ(2*factorial(2*j+2))/QQ(2**(2*j+2))/abs(bern[2*j+2])
for j in range(limit+1)]
t = 1/(pi*pi)
values=[]
for k in range(limit+1):
v=field(0)
for j in reversed(range(k+1)):
v=v*t+field((-1)**j*comb(k,j)*a[j])
values.append(v*t)
return values, pi, {'precision_bits':precision,'arctan_terms_5':n5,
'arctan_terms_239':n239,'bernoulli_max_index':2*limit+2,
'bernoulli_controls':{str(k):str(v) for k,v in controls.items()}}
def assert_narrow(value, bits):
assert value.is_finite(), 'Nonfinite enclosure'
f=value.parent()
assert f(value.rad()) < abs(f(value.mid()))*f(QQ(1)/QQ(2**bits)), \
'Enclosure has insufficient relative accuracy'
def exact_q(decimal):
value=Fraction(decimal)
return QQ(value.numerator)/QQ(value.denominator)
def text_interval(text, precision=1024):
claim=_claim(text)
assert claim is not None
center,radius=map(exact_q,claim)
return RealIntervalField(precision)(center-radius,center+radius),center,radius
def contains_text(text, value):
interval,_,_=text_interval(text)
return interval.lower() <= value.lower() and value.upper() <= interval.upper()
def identities(values, minimum_bits):
field=RealBallField(1024)
values=[field(v) for v in values]
pi=field.pi()
controls=[6/pi**2,6/pi**2-90/pi**4,
6/pi**2-180/pi**4+945/pi**6]
for k,expected in enumerate(controls):
assert_narrow(values[k],minimum_bits)
assert values[k].overlaps(expected)
inversion=[]
egf=[]
f=[1/field(2*j+2).zeta() for j in range(len(values))]
for k in range(len(values)):
inverse=sum(((-1)**j*comb(k,j)*values[j] for j in range(k+1)),field(0))
assert_narrow(inverse,minimum_bits)
assert inverse.overlaps(f[k]), ('binomial inversion',k)
inversion.append(int(inverse.accuracy()))
convolution=sum((field(QQ((-1)**j)/QQ(factorial(j)*factorial(k-j)))*f[j]
for j in range(k+1)),field(0))
coefficient=values[k]/factorial(k)
assert_narrow(coefficient,minimum_bits)
assert_narrow(convolution,minimum_bits)
assert coefficient.overlaps(convolution), ('exponential generating function',k)
egf.append(int(convolution.accuracy()))
return {'closed_form_controls':3,'inversion_orders_checked':len(values),
'inversion_min_relative_accuracy_bits':min(inversion),
'egf_coefficients_checked':len(values),
'egf_rhs_min_relative_accuracy_bits':min(egf)}
def run(stored=False):
start=time.monotonic()
generator=BaezDuarte()
assert list(generator.enumerate()) == [{'k':str(k)} for k in range(MAX_K+1)]
values=list(coefficient_balls(100))
reference,pi,meta=independent(MAX_K)
for value in values:assert_narrow(value,360)
for value in reference:assert_narrow(value,450)
assert not exactness(values)
assert not names_its_rings(__file__.replace('checks.py','generate.py'))
assert all(a.overlaps(b) for a,b in zip(values,reference))
texts=[to_text(value,100) for value in values]
assert all(digits_of(text)==100 for text in texts)
assert all(contains_text(text,ref) for text,ref in zip(texts,reference))
# Controls show both the computation and the comparison can reject errors.
assert not values[1].overlaps(-reference[1])
interval,center,radius=text_interval(texts[0])
bad=RealBallField(1024)(center+10*radius)
assert not contains_text(texts[0],bad)
assert values[0] > 0 and all(value < 0 for value in values[1:])
generated_identities=identities(values,100)
evidence={
'table':'T644','mode':'stored' if stored else 'pre-fill',
'skill_manifest':'agents/runs/20261101T072807Z-build-skills/manifest.json',
'skill_sha256_local':'a2b28d9953264d7b0ca62663785c2952190c5be6d6d3597fd695d9a54d606278',
'skill_sha256_live':'a2b28d9953264d7b0ca62663785c2952190c5be6d6d3597fd695d9a54d606278',
'versions':{'sage':SAGE_VERSION,'numberdb':numberdb.__version__,'python':sys.version.split()[0]},
'entries':len(values),'digits':100,'generator_working_bits':values[0].parent().precision(),
'generator_min_relative_accuracy_bits':min(int(x.accuracy()) for x in values),
'independent_min_relative_accuracy_bits':min(int(x.accuracy()) for x in reference),
'independent_method':meta,'independent_enclosures_overlapping':len(values),
'serialized_intervals_containing_independent_enclosures':len(values),
'controls':['Machin tangent identity exact','B_0,B_1,B_2,B_4,B_6,B_12 exact','c_0,c_1,c_2 closed forms',
'c_1 sign reversal rejected','ten-last-place perturbation rejected'],
'generated_identities':generated_identities,
'raw_measurement':measure(dict(enumerate(values))),
'serialized_measurement':measure(dict(enumerate(texts))),
'values':[{'k':str(k),'number':text} for k,text in enumerate(texts)]}
if stored:
doc=numberdb.table('T644')
entries=_stored_entries(doc)
assert list(entries)==[str(k) for k in range(MAX_K+1)]
stored_values=[]
for k,ref in enumerate(reference):
text=entries[str(k)]
assert digits_of(text)==100
assert contains_text(text,ref), ('stored enclosure',k)
interval,_,_=text_interval(text)
stored_values.append(RealBallField(1024)(interval))
evidence['stored_containment_checks']=len(stored_values)
evidence['stored_identities']=identities(stored_values,60)
report=generator.verify(sample=None)
assert report.ok
evidence['verify']={'checked':report.checked,'matched':report.matched,'ok':report.ok}
evidence['elapsed_seconds']=round(time.monotonic()-start,3)
return evidence
if __name__=='__main__':
configure_key_from_stdin()
print(json.dumps(run('--stored' in sys.argv or os.environ.get('CHECK_STORED')=='1'),indent=2))