back to table · edit · history · where entries came from · files · download
3989 bytes, as of the version from 2026-10-27 16:59 (current). Recorded here, not run.
"""Checks supporting T553 presentation repairs; no writes.
Run: sage -python repair_checks.py, alongside generate.py and
repair-inputs.json.
The archived manifest and its hashes are retained in every report.
"""
import json
from decimal import Decimal
from fractions import Fraction
from pathlib import Path
import importlib.metadata
import numberdb.sage as numberdb
from sage.rings.rational_field import QQ
from sage.rings.real_arb import RealBallField
from sage.rings.complex_arb import ComplexBallField
from sage.env import SAGE_VERSION
from generate import PolygonCapacity
P=Path(__file__).parent
load=lambda name:json.loads((P/name).read_text())
inputs=load('repair-inputs.json')
doc=inputs['table']
R=RealBallField(700)
def rational(s):
f=Fraction(s)
return QQ(f.numerator)/QQ(f.denominator)
def check_decimal(ball,text):
assert ball.is_finite()
centre=rational(text); ulp=QQ(10)**int(Decimal(text).as_tuple().exponent)
lo,hi=QQ(ball.lower()),QQ(ball.upper())
assert 0 < hi-lo < QQ(1)/10**130
assert centre-ulp <= lo <= hi <= centre+ulp
return str(R(max(abs(lo-centre),abs(hi-centre))/ulp))
def value(record): return record['number'] if isinstance(record,dict) else record
# Known answer control, followed by a monotone arithmetic/geometric bracket.
def agm_bracket(x,y):
a,b=max(x,y),min(x,y)
for _ in range(12): a,b=(a+b)/2,(a*b).sqrt()
lo,hi=QQ(b.lower()),QQ(a.upper())
assert 0 <= hi-lo < QQ(1)/10**190
return R((lo+hi)/2,(hi-lo)/2)
assert agm_bracket(R(1),R(1)).contains_exact(1)
agm=agm_bracket(R(1),1/R(2).sqrt())
reverse=agm_bracket(1/R(2).sqrt(),R(1))
assert agm.overlaps(reverse)
agm_text=inputs['agm_table']['number']
agm_ulp=check_decimal(agm,agm_text)
square=1/(R(2).sqrt()*agm)
square_ulp=check_decimal(square,value(doc['Numbers']['4']))
k=ComplexBallField(700)(QQ(1)/2).elliptic_k()
assert k.imag().contains_zero()
k=k.real(); k_ulp=check_decimal(k,inputs['elliptic_table']['number'])
assert (R(2).sqrt()/R.pi()*k).overlaps(square)
g=PolygonCapacity(); parameters=list(g.enumerate())
assert [str(p['n']) for p in parameters]==list(doc['Numbers'])
ns={};exec(doc['Programs']['sage']['code'],ns)
assert list(ns['values'])==list(range(3,31))
rows=[]
for p in parameters:
n=p['n']; v=g.value(p,100)
err=check_decimal(v,value(doc['Numbers'][str(n)]))
assert v.overlaps(ns['values'][n])
a=QQ(1)/QQ(n)
alternate=n*(R.pi()*a).sin()*R(1+a).gamma()/(R(2)**(2*a)*R.pi().sqrt()*R(QQ(1)/2+a).gamma())
assert v.overlaps(alternate)
check_decimal(alternate,value(doc['Numbers'][str(n)]))
area=R(n)/2*(2*R.pi()/n).sin()
area_row=value(inputs['area_table']['numbers'][str(n)])
# Existing areas have fewer digits; require their intervals to contain our ball.
center=rational(area_row);ulp=QQ(10)**int(Decimal(area_row).as_tuple().exponent)
assert center-ulp <= QQ(area.lower()) <= QQ(area.upper()) <= center+ulp
assert (R.pi()/n).cos() < v < 1
rows.append({'n':n,'decimal_error_ulps_upper':err,'gamma_identity':'pass','area_normalization':'pass'})
report={'table':'T553','scope':'Presentation repair verification; no new or changed values.',
'policy':inputs['policy'],
'sources':[{'url':'https://dlmf.nist.gov/19.8#i','locators':'19.8.1, symmetry paragraph, 19.8.4 and 19.8.5; modulus converted to parameter m=k^2'},
{'url':'https://dlmf.nist.gov/19.2#E8','locators':'19.2.4 and 19.2.8'}],
'argument_uncertainty':'Exact integer/rational inputs; Arb encloses sqrt(2), Gamma arguments and subsequent arithmetic throughout.',
'control':'AGM(1,1)=1 passed before the square check.',
'agm_t51_error_ulps_upper':agm_ulp,'square_capacity_error_ulps_upper':square_ulp,
'elliptic_k_t25_error_ulps_upper':k_ulp,'agm_enclosure_width':str(R(QQ(agm.upper())-QQ(agm.lower()))),
'generator_and_program_enumeration':'Both reproduce exactly n=3,...,30.', 'rows':rows,
'versions':{'sage':SAGE_VERSION,'numberdb':importlib.metadata.version('numberdb')}}
print(json.dumps(report,indent=2))