back to table · edit · history · where entries came from · files · download
4507 bytes, as of the version from 2026-10-29 02:59 (current). Recorded here, not run.
"""Internal solid angles of regular simplices -- numberdb.org/T598.
Run with SageMath (install numberdb with sage -pip install numberdb):
sage -python generate.py # verify every stored entry
sage -python generate.py --fill # fill this draft, NUMBERDB_API_KEY required
sage -python checks.py # independent local numerical checks
Filling entries does not publish or offer the draft. The default enumeration
contains all 77 pairs 2 <= d <= 12, 0 <= k < d. See argument-selection.md for
its justification and numerical-method.md for the proof of the tail bound.
The first two rational cases are decided geometrically, never from digits.
The remaining values use validated Arb quadrature of the entire, symmetrized
imaginary-normal integrand. All arithmetic, including arguments, is exact or
ball arithmetic. WORKING_GUARD is 167 bits; the separate quadrature target is
97 bits beyond numberdb.bits(digits). The complete run measures retention.
"""
import os
import sys
from functools import lru_cache
from math import isqrt
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
WORKING_GUARD = 167
QUADRATURE_GUARD = 97
MAX_DIMENSION = 12
@lru_cache(maxsize=None)
def quadrature(d, k, digits=100):
"""Return a certified angle, even for the exact control cases."""
if not (d >= 2 and 0 <= k < d):
raise ValueError('Require d >= 2 and 0 <= k < d')
m = d - k
n = d + 1
bits = numberdb.bits(digits, losing=WORKING_GUARD)
C = ComplexBallField(bits)
R = RealBallField(bits)
scale = C(0, 1) / C(2 * n).sqrt()
normalizer = C(2 * C.pi()).sqrt()
# a*T^2/2 >= 3*(digits+30), with exact integer arithmetic.
cutoff = isqrt((6 * n * (digits + 30) + k) // (k + 1)) + 1
tolerance = R(2) ** (-numberdb.bits(digits, losing=QUADRATURE_GUARD))
def integrand(x, analytic):
# Entire in x. Real-part extraction here would break analyticity.
e = (scale * x).erf()
return ((1 + e) ** m + (1 - e) ** m) * (-x * x / 2).exp() / (2 ** m * normalizer)
integral = C.integral(integrand, 0, cutoff,
abs_tol=tolerance, rel_tol=tolerance)
assert abs(integral.real()) < 1 and integral.imag().contains_zero()
assert abs(integral.imag()) < R(10) ** (-digits)
a = R(QQ(k + 1) / QQ(n))
tail = (R(2) / R.pi()).sqrt() * (-a * cutoff * cutoff / 2).exp() / (a * cutoff)
result = integral.real().add_error(tail.upper())
assert result > 0 and result < 1
assert result.rad_as_ball() / abs(result) < R(10) ** (-(digits + 10))
return result
class InternalSimplexAngles(numberdb.Generator):
table = 'T598'
parameters = ('d', 'k')
type = 'R'
digits = 100
rigour = 'proven'
files = ('generate.py', 'checks.py', 'check.py', 'argument-selection.md',
'numerical-method.md', 'numerical-evidence.json', 'program.py')
def enumerate(self, max_dimension=MAX_DIMENSION):
for d in range(2, max_dimension + 1):
for k in range(d):
yield {'d': d, 'k': k}
def value(self, params, digits):
d, k = int(params['d']), int(params['k'])
if not (d >= 2 and 0 <= k < d):
raise ValueError('Require d >= 2 and 0 <= k < d')
if k == d - 1:
return QQ(1) / QQ(2)
if (d, k) == (2, 0):
return {'number': QQ(1) / QQ(6), 'comment': 'A vertex of an equilateral triangle has angle $\\pi/3$, which is $1/6$ of the full circle.'}
value = quadrature(d, k, digits)
if (d, k) == (3, 0):
return {'number': value, 'comment': 'The solid angle at a vertex of a regular tetrahedron is $4\\pi\\beta_{3,0}=3\\arccos(1/3)-\\pi$ steradians.'}
return value
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
generator = InternalSimplexAngles()
if '--fill' in sys.argv or os.environ.get('NUMBERDB_PUBLISH') == '1':
print(generator.publish(message='Fill 77 certified internal simplex angles, with independent geometric checks.'))
else:
report = generator.verify(sample=None)
print(report)
sys.exit(0 if report.ok else 1)