back to table · edit · history · where entries came from · files · download
8903 bytes, as of the version from 2026-10-29 02:59 (current). Recorded here, not run.
"""Independent checks for internal simplex angles -- numberdb.org/T598.
sage -python checks.py
sage -python checks.py --stored # authenticated key needed while a draft
The main cross-check uses real Gaussian quadrature for external angles and
conic Gauss-Bonnet to recover internal angles. It does not use the generator's
imaginary-normal integrand. Closed forms, exact projected Gram matrices and
Gram's angle-sum relation supply further geometric checks. Full derivations
and the archived skill manifest are recorded in numerical-method.md.
"""
import json
import os
import sys
import time
from functools import lru_cache
from fractions import Fraction
from math import comb, isqrt
from pathlib import Path
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 numberdb._write import to_text
from numberdb._generate import _stored_entries
from numberdb._compare import _claim
from check import exactness, measure, agrees_with, names_its_rings
from generate import InternalSimplexAngles, quadrature
CHECK_BITS = 630
C = ComplexBallField(CHECK_BITS)
R = RealBallField(CHECK_BITS)
def sharp(value, digits=100):
value = R(value)
assert abs(value) < 2, 'Nonfinite or implausible value'
assert value.rad_as_ball() < R(10) ** (-digits) * abs(value), 'Imprecise check'
return value
def overlaps(a, b, digits=100):
return sharp(a, digits).overlaps(sharp(b, digits))
@lru_cache(maxsize=None)
def external_quadrature(m, q):
"""g_m(1/q), via a nonnegative real-normal integrand at 630 bits."""
scale = C(2 * q).sqrt()
normalization = C(2 * C.pi()).sqrt()
tolerance = R(2) ** (-560)
cutoff = 34
def f(x, analytic):
e = (x / scale).erf()
return ((1 + e) ** m + (1 - e) ** m) * (-x * x / 2).exp() / (2 ** m * normalization)
z = C.integral(f, 0, cutoff, abs_tol=tolerance, rel_tol=tolerance)
assert z.imag().contains_zero() and abs(z.imag()) < tolerance
# For real x each Phi power is at most 1; both Gaussian tails.
tail = 2 * (-R(cutoff) ** 2 / 2).exp() / (R(cutoff) * (2 * R.pi()).sqrt())
return sharp(z.real().add_error(tail.upper()), 150)
@lru_cache(maxsize=None)
def external(d, k):
if d == k:
return QQ(1)
if k == 0:
return QQ(1) / QQ(d + 1)
if k == d - 1:
return QQ(1) / QQ(2)
return external_quadrature(d - k, k + 1)
@lru_cache(maxsize=None)
def independent(d, k):
"""Recover beta using one parity of conic intrinsic volumes."""
if d == k:
return R(1)
if d == k + 1:
return R(QQ(1) / QQ(2))
m = d - k
answer = R(QQ(1) / QQ(2))
for j in range(m % 2, m, 2):
answer -= comb(m, j) * independent(k + j, k) * external(d, k + j)
return sharp(answer, 140)
def stored_values():
doc = numberdb.table('T598')
written = _stored_entries(doc)
values = {}
for key, text in written.items():
d, k = map(int, key.split(','))
if '/' in text or ('.' not in text and 'e' not in text.lower()):
values[d, k] = QQ(text)
else:
center, radius = _claim(text)
c, r = Fraction(center), Fraction(radius)
values[d, k] = R(QQ(c.numerator) / c.denominator).add_error(
R(QQ(r.numerator) / r.denominator).upper())
return values, written
def run(stored=False):
start = time.monotonic()
# Known nonzero controls run before any family comparisons.
for m, q, known in [(0, 1, QQ(1)), (1, 1, QQ(1)/2),
(2, 1, QQ(1)/3), (3, 1, QQ(1)/4)]:
assert overlaps(external_quadrature(m, q), known, 145)
assert overlaps(quadrature(2, 0), QQ(1)/6, 110)
assert overlaps(quadrature(2, 1), QQ(1)/2, 110)
print('PASS: 6 nonzero quadrature controls', flush=True)
gen = InternalSimplexAngles()
keys = [(p['d'], p['k']) for p in gen.enumerate()]
assert keys == [(d, k) for d in range(2, 13) for k in range(d)]
assert len(keys) == 77
produced = {}
written = {}
for d, k in keys:
value = gen.value({'d':d, 'k':k}, 100)
value = value['number'] if isinstance(value, dict) else value
produced[d, k] = value
written['%d,%d' % (d,k)] = to_text(value, digits=100)
assert not exactness(produced)
assert not names_its_rings(str(Path(__file__).with_name('generate.py')))
relative_radii = {str(key):float(R(v.rad_as_ball()) / abs(R(v)))
for key,v in produced.items() if hasattr(v,'rad_as_ball')}
print('PASS: enumeration, exactness and named imports', flush=True)
for d, k in keys:
high = quadrature(d, k, 140)
assert overlaps(produced[d,k], high, 110)
print('PASS: all 77 entries agree with 140-digit recomputation', flush=True)
for d, k in keys:
expected = independent(d,k)
assert overlaps(produced[d,k], expected, 110), (d,k)
assert not agrees_with(produced, lambda key: independent(*key))
print('PASS: all 77 entries agree with external-angle/Gauss-Bonnet computation', flush=True)
values = stored_values()[0] if stored else produced
assert list(values) == keys
# Written intervals use one unit in their last place, so 97 is the
# diagnostic width threshold for values displayed with 100 digits.
threshold = 97 if stored else 110
count = {'facets':0, 'codimension_two':0, 'codimension_three':0,
'exact_gram_matrices':0, 'gram_angle_sums':0,
'gauss_bonnet_recurrences':0}
for d,k in keys:
v = sharp(values[d,k], threshold)
assert v > 0 and v < 1
assert overlaps(v, independent(d,k), threshold)
m=d-k
# Project e_j - e_0 orthogonally off the face direction space:
# u_j = e_j - mean(e_0,...,e_k). Build their actual coordinates.
vectors = [[(QQ(1) if h==j else QQ(0)) -
(QQ(1)/QQ(k+1) if h<=k else QQ(0))
for h in range(d+1)] for j in range(k+1,d+1)]
gram = [[sum(x*y for x,y in zip(u,w)) for w in vectors] for u in vectors]
inv = [[QQ(int(i==j))-QQ(1)/QQ(d+1) for j in range(m)] for i in range(m)]
for i in range(m):
for j in range(m):
assert gram[i][j] == QQ(int(i==j))+QQ(1)/QQ(k+1)
assert sum(gram[i][h]*inv[h][j] for h in range(m)) == int(i==j)
count['exact_gram_matrices']+=1
if m==1:
assert overlaps(v, QQ(1)/2, threshold); count['facets']+=1
if m==2:
assert overlaps(v, (R(1)/d).arccos()/(2*R.pi()), threshold)
count['codimension_two']+=1
if m==3:
assert overlaps(v, QQ(1)/8 - 3*(R(1)/d).arcsin()/(4*R.pi()), threshold)
count['codimension_three']+=1
# Evaluate the full parity identity using the chosen values themselves.
total=R(v)
for j in range(m%2,m,2):
inner=R(1) if j==0 else (R(QQ(1)/2) if j==1 else R(values[k+j,k]))
total += comb(m,j)*inner*external(d,k+j)
assert (total-QQ(1)/2).contains_zero() and total.rad_as_ball()<R(10)**(-(threshold-8))
count['gauss_bonnet_recurrences']+=1
for d in range(2,13):
total=sum((-1)**k*comb(d+1,k+1)*R(values[d,k]) for k in range(d))
assert (total-(-1)**(d+1)).contains_zero()
assert total.rad_as_ball()<R(10)**(-(threshold-8))
count['gram_angle_sums']+=1
assert overlaps(values[2,0], QQ(1)/6, threshold)
assert overlaps(4*R.pi()*R(values[3,0]), 3*(R(1)/3).arccos()-R.pi(), threshold)
report={'table':'T598','mode':'stored' if stored else 'local','entries':77,
'exact_rationals':12,'approximate_values':65,'digits':100,
'working_bits':numberdb.bits(100,losing=167),
'independent_bits':CHECK_BITS,'higher_precision_digits':140,
'checks':count,'max_relative_radius':max(relative_radii.values()),
'worst_relative_radius_at':max(relative_radii,key=relative_radii.get),
'measurement':measure(written),'environment':gen.environment(),
'skill_manifest':'agents/runs/20261029T024405Z-build-skills/manifest.json',
'elapsed_seconds':time.monotonic()-start}
print('EVIDENCE_JSON '+json.dumps(report,sort_keys=True),flush=True)
print('VALUES_JSON '+json.dumps(written,sort_keys=False),flush=True)
print('PASS: all geometric identities, covariance and precision checks',flush=True)
return report
if __name__=='__main__':
if os.environ.get('NUMBERDB_KEY_FROM_STDIN')=='1':
key=sys.stdin.read().strip()
if '=' in key and key.split('=',1)[0].replace('_','').isupper():
key=key.split('=',1)[1].strip().strip('\"\'')
os.environ['NUMBERDB_API_KEY']=key
run(stored='--stored' in sys.argv or os.environ.get('CHECK_STORED')=='1')