check_example.py

back to table · edit · history · where entries came from · files · download

3325 bytes, as of the version from 2026-10-30 04:30 (current). Recorded here, not run.

"""Check T617's reader program under Sage; no website writes.

Run: agents/sage.sh check_example.py pressure_example.py generate.py solver.py
     radial_check.py T617-before.json skill-evidence.json
Independent references for new examples use the existing radial pressure
integral calculation, not the SciPy surface-derivative formula.
"""
import json, sys, hashlib, time
from pathlib import Path
import numberdb.sage as numberdb
from sage.rings.real_mpfr import RealField
from sage.rings.rational_field import QQ
from sage.env import SAGE_VERSION
import scipy
from pressure_example import pressure
from generate import LaneEmdenPressure
from radial_check import radial
root=Path(__file__).parent;F=RealField(240)
doc=json.loads((root/'T617-before.json').read_text())
print(json.dumps({'skill':json.loads((root/'skill-evidence.json').read_text()),'sage':SAGE_VERSION,'scipy':scipy.__version__,'python':sys.version,'input_sha256':hashlib.sha256((root/'T617-before.json').read_bytes()).hexdigest(),'program_sha256':hashlib.sha256((root/'pressure_example.py').read_bytes()).hexdigest()}),flush=True)
g=LaneEmdenPressure()
for n in ['13/4','1/3']:
    try:g.value({'n':n},100)
    except ValueError as e:print(json.dumps({'original_program_rejection':n,'message':str(e)}),flush=True)
    else:raise AssertionError('Original failure no longer reproduces')
assert [p['n'] for p in g.enumerate()]==list(doc['Numbers'])
print(json.dumps({'enumeration_matches':True,'count':len(doc['Numbers'])}),flush=True)
worst=F(0)
for n,entry in doc['Numbers'].items():
    value=entry['number'] if isinstance(entry,dict) else entry
    ref=F(value); actual=F(pressure(n));tight=F(pressure(n,rtol=3e-14,atol=3e-16,start=5e-5))
    rel=abs(actual-ref)/abs(ref);convergence=abs(actual-tight)/abs(tight)
    assert rel<F('1e-8') and convergence<F('1e-8')
    worst=max(worst,rel)
    print(json.dumps({'n':n,'reference':'stored 100-digit value','actual':str(actual),'relative_error':str(rel),'tolerance_start_change':str(convergence)}),flush=True)
for n,ref in [('0',3/(8*F.pi())),('1',F.pi()/8)]:
    assert abs(F(pressure(n))/ref-1)<F('1e-8')
    assert abs(F(pressure(n))*(1+F('1e-5'))/ref-1)>F('1e-8')
print(json.dumps({'analytic_controls_pass':True,'displaced_controls_rejected':True}),flush=True)
for n in ['13/4','1/3']:
    started=time.monotonic();mass,meta=radial(n)
    ref=F(meta['radius'])**3/(4*F.pi()*(5-F(QQ(n)))*F(meta['pressure_integral']))
    a=F(pressure(n));b=F(pressure(n,rtol=3e-14,atol=3e-16,start=5e-5))
    rel=abs(a/ref-1);assert rel<F('1e-8') and abs(b/ref-1)<F('1e-8')
    print(json.dumps({'n':n,'reference':'70-digit radial Taylor pressure integral; heuristic endpoint correction','reference_value':str(ref),'relative_error':str(rel),'tight_relative_error':str(abs(b/ref-1)),'metadata':meta,'seconds':time.monotonic()-started}),flush=True)
for n in ['-1','5','6']:
    try:pressure(n)
    except ValueError:pass
    else:raise AssertionError('Domain validation failed')
try:pressure('3',limit=1)
except RuntimeError:pass
else:raise AssertionError('Missing surface accepted')
print(json.dumps({'ok':True,'stored_values_changed':False,'stored_worst_relative_error':str(worst),'accuracy_claim':'8 displayed significant digits, empirically checked, not certified; no table values computed by this example'}),flush=True)