back to table · edit · history · where entries came from · files · download
4341 bytes, as of the version from 2026-10-30 02:35 (current). Recorded here, not run.
"""Check stored T615 intervals, analytic anchors and companion identities.
Run: sage -python readback_checks.py, with NUMBERDB_API_KEY for draft reads.
The campaign wrapper supports NUMBERDB_KEY_FROM_STDIN=1. Companion tables are
read through the public API with this account's draft authorization. The
repository check.stored helper requires Django and is not used on this host.
"""
import os,sys,json,pathlib
from decimal import Decimal,localcontext
import numberdb.sage as numberdb
from sage.rings.real_arb import RealBallField
from generate import INDICES
def number(item):
return str(item['number'] if isinstance(item,dict) else item)
def interval(text):
d=Decimal(text)
if '.' not in text and 'e' not in text.lower():return d,d,Decimal(0)
ulp=Decimal(10)**d.as_tuple().exponent
return d-ulp,d+ulp,ulp
def main():
if os.environ.get('NUMBERDB_KEY_FROM_STDIN')=='1':
token=sys.stdin.read().strip()
if '=' in token and token.split('=',1)[0].replace('_','').isalnum():
token=token.split('=',1)[1].strip().strip('\"\'')
os.environ['NUMBERDB_API_KEY']=token
d=numberdb.table('T615');radii=numberdb.table('T613');masses=numberdb.table('T614')
assert d['Title']=='Central-to-mean density ratio of Lane–Emden polytropes'
assert d['Display properties']['number-header']==r'$C_n=\rho_c/\bar\rho$'
assert d['Parameters']['n']['type']=='Q'
assert list(d['Numbers'])==list(INDICES)
p=pathlib.Path(__file__).parent
computed=json.loads((p/'numerical-evidence.json').read_text())['values']
independent=json.loads((p/'independent-evidence.json').read_text())['values']
skill=json.loads((p/'skill-evidence.json').read_text())
assert skill['match'] is True
assert set(skill['versions'].values())=={'28ecbeeccb0c5e5d6eb541b5a5b0f88aec47b154374d909085116d34f5629199'}
report={'table':'T615','skill_manifest':skill['manifest'],'rows':{}}
B=RealBallField(700)
with localcontext() as ctx:
ctx.prec=250
for n in INDICES:
s=number(d['Numbers'][n]);lo,hi,ulp=interval(s);center=Decimal(s)
if n=='0':
assert s=='1'
else:
assert len(center.as_tuple().digits)>=100,(n,'too few digits')
if n=='1':
target=B.pi()**2/3
assert target.is_finite() and B(target.rad())<B('1e-190')
assert B(str(lo))<target and target<B(str(hi))
for work in ['150','210']:
assert lo<Decimal(computed[n][work]['ratio'])<hi,(n,work)
error=abs(Decimal(independent[n]['ratio'])-center)/abs(center)
assert error<Decimal('1e-40'),(n,error)
rlo,rhi,_=interval(number(radii['Numbers'][n]))
mlo,mhi,_=interval(number(masses['Numbers'][n]))
assert rlo<Decimal(computed[n]['210']['radius'])<rhi,('radius',n)
assert mlo<Decimal(computed[n]['210']['mass'])<mhi,('mass',n)
# Positive interval arithmetic: the stored companion intervals must
# give a concentration interval intersecting this row's interval.
clo=rlo**3/(3*mhi);chi=rhi**3/(3*mlo)
assert max(clo,lo)<=min(chi,hi),('companion identity',n)
# A separate integral identity using the stored concentration to
# infer mass from the independently integrated radius.
xi=Decimal(independent[n]['radius']);J=Decimal(independent[n]['pressure_integral'])
mu=xi**3/(3*center)
nn=Decimal(n.split('/')[0])/(Decimal(n.split('/')[1]) if '/' in n else 1)
virial=(nn+1)*mu*mu/((5-nn)*xi)
verr=abs(J-virial)/abs(virial)
assert verr<Decimal('1e-40'),('stored virial',n,verr)
report['rows'][n]={'exact':n=='0','stored_interval_pass':True,
'significant_digits':None if n=='0' else len(center.as_tuple().digits),
'independent_relative_difference':str(error),'companion_identity_pass':True,
'stored_virial_relative_error':str(verr)}
anchor=number(d['Numbers']['1']);lo,hi,ulp=interval(anchor)
assert not lo<Decimal(anchor)+100*ulp<hi
report['perturbed_value_rejected']=True
print('DATA '+json.dumps(report),flush=True)
if __name__=='__main__':main()