checks.py

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

5593 bytes, as of the version from 2026-10-25 08:29 (current). Recorded here, not run.

"""Independent numerical checks for Tammes separations -- numberdb.org/T479.

Run: sage -python checks.py > checks-result.json
Requires the attached generate.py, check.py and sources.json beside this file.
Optional stored-document.json makes the same checks on API readback values.
Coordinates are radially normalized. This verifies separation, not global
optimality, and does not rerun the source's search for spherical codes.
"""
import json
from pathlib import Path
from decimal import Decimal
import numberdb.sage as numberdb
from sage.rings.real_arb import RealBallField
from sage.rings.rational_field import QQ
from numberdb._write import to_text
from check import exactness, measure, names_its_rings
import generate

R=RealBallField(768)

def rat(text):
    a,b=Decimal(str(text)).as_integer_ratio()
    return QQ(a)/QQ(b)

def ball(text):
    if '.' not in text and 'e' not in text.lower():return R(QQ(text))
    ulp=QQ(10)**Decimal(text).as_tuple().exponent
    return R(rat(text),ulp)

def finite(x):
    assert x.is_finite()
    assert R(x.rad())<R('1e-7'),(str(x),str(x.rad()))

def separation(rows):
    points=[]
    normerr=R(0)
    for raw in rows:
        v=[R(rat(x)) for x in raw]
        norm=sum(x*x for x in v)
        err=abs(norm-1)
        if err.mid()>normerr.mid():normerr=err
        points.append([x/norm.sqrt() for x in v])
    dots=[];distances=[]
    for i,u in enumerate(points):
        for v in points[:i]:
            dots.append(sum(a*b for a,b in zip(u,v)))
            distances.append(sum((a-b)**2 for a,b in zip(u,v)).sqrt())
    c=max(dots,key=lambda v:v.mid())
    d=min(distances,key=lambda v:v.mid())
    finite(c);finite(d)
    assert abs(d*d-(2-2*c))<R('1e-60')
    return c,d,normerr

# Controls first: the computation must recover known nonzero separations.
controls=[([['1','0','0'],['-1','0','0']],R(-1),R(2)),
          ([['1','1','1'],['1','-1','-1'],['-1','1','-1'],['-1','-1','1']],-R(1)/3,(R(8)/3).sqrt()),
          ([['1','0','0'],['-1','0','0'],['0','1','0'],['0','-1','0'],['0','0','1'],['0','0','-1']],R(0),R(2).sqrt())]
for points,c0,d0 in controls:
    c,d,_=separation(points)
    assert abs(c-c0)<R('1e-60') and abs(d-d0)<R('1e-60')

G=generate.Tammes()
values={};entries={};evidence=[]
for params in G.enumerate():
    key=f"{params['n']},{params['normalisation']}"
    v=G.value(params,100); values[key]=v
    text=to_text(v['number'],digits=v['digits'])
    row={'number':text}
    if v.get('comment'):row['comment']=v['comment']
    entries.setdefault(str(params['n']),{})[params['normalisation']]=row
assert not exactness(values)
assert not names_its_rings(str(Path(generate.__file__)))
readback=Path(__file__).with_name('stored-document.json')
if readback.exists():
    stored=json.loads(readback.read_text())['Numbers']
    assert set(stored)==set(entries)
    for n in entries:
        assert set(stored[n])==set(entries[n])
        for norm in entries[n]:
            assert stored[n][norm]['number']==entries[n][norm]['number'],(n,norm)
    entries=stored

minalgebraic=None
for n in range(4,101):
    row=generate.DATA['rows'][str(n)]
    c,d,normerr=separation(generate.DATA['coordinates'][str(n)]['points'])
    expected=R(generate.cosine(n))
    tolerance=R('1e-8') if n<=6 else R('1e-12') if 'polynomial' in row else R('2e-12')
    assert abs(c-expected)<tolerance,(n,'coordinate cosine',str(c-expected))
    if n>=7:
        assert abs(c-R(rat(row['cosine'])))<R('1e-12'),(n,'source cosine',str(c))
    a=ball(entries[str(n)]['angle']['number'])
    b=ball(entries[str(n)]['chord']['number'])
    finite(a);finite(b)
    # The two normalisations must agree as stored intervals, in both directions.
    assert b.overlaps(2*(a*R.pi()/360).sin()),(n,'angle to chord')
    assert a.overlaps(360*(b/2).arcsin()/R.pi()),(n,'chord to angle')
    assert abs(a-c.arccos()*180/R.pi())< (R('1e-6') if n<=6 else R('1e-7'))
    assert abs(b-d)< (R('1e-8') if n<=6 else R('1e-9'))
    if 'polynomial' not in row and n>=7:
        assert a.overlaps(c.arccos()*180/R.pi()),(n,'coordinate angle outside stored interval')
        assert b.overlaps(d),(n,'coordinate chord outside stored interval')
    if n<=6 or 'polynomial' in row:
        for norm in ['angle','chord']:
            z=values[f'{n},{norm}']['number']
            if not isinstance(z,int):
                acc=int((-R(z.rad()/abs(z.mid())).log()/R(2).log()).lower())
                minalgebraic=acc if minalgebraic is None else min(minalgebraic,acc)
                assert acc>332
    # Sloane is an older independent tabulation: differences are recorded.
    da=abs(c.arccos()*180/R.pi()-R(rat(row['sloane_angle'])))
    evidence.append({'n':n,'coordinate_minus_expected_cosine':str(c-expected),'squared_norm_error':str(normerr),'difference_from_sloane_degrees':str(da),'sloane_agrees_7_decimals':bool(da<R('1e-7'))})

# Outside the optimized-coordinate family: known regular-polyhedron formulas.
assert ball(entries['4']['chord']['number']).overlaps((R(8)/3).sqrt())
assert entries['5']['angle']['number']=='90' and entries['6']['angle']['number']=='90'
assert ball(entries['12']['chord']['number']).overlaps((2-2/R(5).sqrt()).sqrt())
assert ball(entries['8']['chord']['number']).overlaps((2-2*(2*R(2).sqrt()-1)/7).sqrt())
assert ball(entries['9']['chord']['number']).overlaps(2/R(3).sqrt())

result={'ok':True,'readback':readback.exists(),'controls':3,'configurations':97,'entries_checked':194,'algebraic_min_relative_accuracy_bits':minalgebraic,'measurement':measure(values),'coordinate_checks':evidence,'Numbers':entries}
print(json.dumps(result,indent=2))