back to table · edit · history · where entries came from · files · download
7065 bytes, as of the version from 2026-10-26 20:13 (current). Recorded here, not run.
"""Independent numerical and stored-interval checks for numberdb.org/T526.
Run: sage -python checks.py [--stored]
The stored check needs NUMBERDB_API_KEY for the private draft.
This file derives polynomial killed-motion resolvents independently of the
K/Q generator; see derivation.md. check.py is the campaign's unmodified helper.
"""
import os
import sys
import json
from fractions import Fraction as F
from math import factorial
from pathlib import Path
import numberdb.sage as numberdb
from sage.rings.rational_field import QQ
from sage.rings.real_arb import RealBallField
from numberdb._write import to_text
from check import exactness, measure, agrees_with, names_its_rings
from generate import MeanderMoments, auxiliary, rational_factor, WORKING_GUARD
def trim(p):
p = list(p)
while len(p)>1 and p[-1]==0: p.pop()
return p
def deriv(p):
return trim([F(i)*p[i] for i in range(1,len(p))] or [F(0)])
def add(p,q,scale=F(1)):
r=[F(0)]*max(len(p),len(q))
for i,x in enumerate(p): r[i]+=x
for i,x in enumerate(q): r[i]+=scale*x
return trim(r)
def poly_coeff(p,i):
return p[i] if i<len(p) else F(0)
def resolvent_slopes(last):
p,r=[F(1)],[F(-1)]
result=[F(1)]
for n in range(1,last+1):
source_p=[F(0)]+p
source_r=[F(0)]+r
p=[F(0)]; term=source_p
while term != [0]:
p=add(p,term); term=deriv(deriv(term))
slope=[F(0)]; term=source_r; divisor=F(2)
while term != [0]:
slope=add(slope,term,1/divisor)
term=deriv(term); divisor*=2
r=[-p[0]]+[slope[i]/F(i+1) for i in range(len(slope))]
r=trim(r)
assert add(p,deriv(deriv(p)),F(-1))==trim(source_p)
assert add([2*x for x in deriv(r)],deriv(deriv(r)),F(-1))==trim(source_r)
assert p[0]+r[0]==0
result.append(poly_coeff(p,1)+poly_coeff(r,1)-r[0])
return result
def independent_factor(n,slope):
"""Cancel Gamma using products, independently of generator factorial ratios."""
base=F(factorial(n))*slope/F(2**(n//2))
if n%2:
for k in range(1,(3*n-1)//2+1): base/=k
else:
for k in range(3*n//2): base/=F(2*k+1,2)
return base
def direct_moment(n,slope,field):
return (field.pi().sqrt()*field(2)**(-field(n)/2)*field(factorial(n))
*field(QQ(slope))/field(QQ(3*n+1)/2).gamma())
def endpoints(text):
text=str(text)
if '.' not in text and 'e' not in text.lower():
q=F(text); return q,q
mant,sep,exp=text.lower().partition('e')
places=len(mant.split('.')[1]) if '.' in mant else 0
ulp=F(10)**((int(exp) if sep else 0)-places)
center=F(text)
return center-ulp,center+ulp
def strict_enclosure(value,lo,hi,field):
assert value.is_finite() and value.accuracy()>=700
assert value > field(QQ(lo)) and value < field(QQ(hi))
def run(stored=False):
assert endpoints('1.25')==(F(124,100),F(126,100))
assert endpoints('12e2')==(F(1100),F(1300))
g=MeanderMoments(); last=12
assert list(g.enumerate())==[{'n':str(n)} for n in range(1,13)]
slopes=resolvent_slopes(last)
rb=RealBallField(768)
assert slopes[0]==1 and slopes[1]==F(3,4)
assert (direct_moment(0,slopes[0],rb)-rb(1)).contains_zero()
controls={0:F(1),1:F(3,4),2:F(59,60),3:F(465,512),4:F(5345,3696)}
for n,want in controls.items():
assert independent_factor(n,slopes[n])==want
factors={n:rational_factor(n) for n in range(last+1)}
assert not exactness(factors)
assert not agrees_with(factors,lambda n:QQ(independent_factor(n,slopes[n])))
k,q=auxiliary(last)
assert not agrees_with(dict(enumerate(q)),lambda n:QQ(slopes[n]))
assert not names_its_rings(str(Path(__file__).with_name('generate.py')))
values={n:g.value({'n':str(n)},100) for n in range(1,last+1)}
assert not exactness(values)
texts={str(n):to_text(v,digits=100) for n,v in values.items()}
if stored:
doc=numberdb.table('T526')
assert list(doc['Numbers'])==list(texts)
texts={n:(v['number'] if isinstance(v,dict) else v) for n,v in doc['Numbers'].items()}
report={'table':'T526','mode':'stored' if stored else 'offline',
'skill_manifest':'agents/runs/20261026T200330Z-build-skills/manifest.json',
'skill_sha256':'b4d51ded0d52cfe3dcc0859de0845c3879b0c2363273ba7e809a7e06ded40b18',
'numberdb_version':numberdb.__version__,
'generator_working_bits':numberdb.bits(100,losing=WORKING_GUARD),
'independent_working_bits':768,
'min_odd_relative_accuracy_bits':min(int(v.accuracy()) for n,v in values.items() if n%2),
'measure':measure(texts), 'controls':'M0 and Table 4 orders 1 through 4; exact interval parser',
'polynomial_ode_residuals':'zero, both equations, all orders 1 through 12',
'coefficient_agreement':'exact, all orders 0 through 12', 'rows':[]}
balls={0:rb(1)}; recovered={0:rb(1)}
for n in range(1,last+1):
text=texts[str(n)]; lo,hi=endpoints(text)
independent=direct_moment(n,slopes[n],rb)
assert independent.is_finite() and independent.accuracy()>=700
if n%2:
strict_enclosure(independent,lo,hi,rb)
# Decimal center written to exactly 100 significant digits.
from decimal import Decimal
assert len(Decimal(text).as_tuple().digits)==100
mid=(lo+hi)/2; rad=(hi-lo)/2
ball=rb(QQ(mid)).add_error(rb(QQ(rad)))
else:
assert lo==hi==independent_factor(n,slopes[n])
assert (independent-rb(QQ(lo))).contains_zero()
ball=rb(QQ(lo))
assert ball.is_finite() and ball.accuracy()>=320
balls[n]=ball
recovered[n]=(ball*rb(QQ(3*n+1)/2).gamma()*rb(2)**(rb(n)/2)
/(rb.pi().sqrt()*rb(factorial(n))))
residual=recovered[n]-rb(QQ(3*n-2)/2)*recovered[n-1]
for j in range(1,n+1): residual-=2*rb(k[j])*recovered[n-j]
assert residual.is_finite() and residual.contains_zero()
assert residual.rad_as_ball()<rb(10)**(-85)
report['rows'].append({'n':n,'Q_n':str(slopes[n]),
'rational_factor':str(independent_factor(n,slopes[n])),
'number':text,'independent_ball_accuracy_bits':int(independent.accuracy()),
'stored_recurrence_residual_radius':str(residual.rad())})
for n in range(1,last):
assert balls[n-1]*balls[n+1]>balls[n]**2
report['log_convexity']='strict for n=1,...,11, including M0=1'
report['all_written_intervals_checked']=True
if stored:
verification=g.verify(sample=None)
assert verification.ok, str(verification)
report['verify']=str(verification)
print(json.dumps(report,indent=2))
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':
numberdb.configure(api_key=sys.stdin.read().strip())
run(stored=stored)