back to table · edit · history · where entries came from · files · download
4608 bytes, as of the version from 2026-10-27 02:01 (current). Recorded here, not run.
"""Independent gamma/elliptic checks for T521's five even additions.
python even-checks.py even-computation.jsonl
No generator or Bessel evaluator is imported. BBBG, arXiv:0801.0891v2,
equations (47), (49), (78), (89), (95), and (100) give the identities used.
The quarter-interval reduction immediately following (81) removes the
logarithmic endpoint from c40. Combining (89) and (95) gives
c42=(-2*c40+3*B1+3*pi^2-3*H)/16. Elliptic K takes the parameter m here,
whereas the paper writes its modulus. Carlson RF avoids subtracting tiny
complementary parameters in the s40 integral. Independent results are
agreement-checked at 140 and 180 digits, not labelled certified enclosures.
"""
import json
import sys
from decimal import Decimal
from pathlib import Path
import mpmath as mp
EVEN_KEYS = ((0, 3, 0), (0, 3, 2), (0, 4, 0), (0, 4, 2), (1, 3, 0))
def reference_values(working):
with mp.workdps(working):
third_gamma = mp.gamma(mp.mpf(1) / 3)
third_zero = 3 * third_gamma ** 6 / (32 * mp.pi * mp.root(4, 3))
third_second = third_zero / 9 - mp.pi ** 4 / (24 * third_zero)
def elliptic_square(angle):
return mp.ellipk(mp.sin(angle) ** 2) ** 2
fourth_zero = 4 * mp.pi * mp.quad(
elliptic_square, [0, mp.pi / 8, mp.pi / 4], method='gauss-legendre')
def subtracted(angle):
if not angle:
return mp.pi ** 3 / 2
return (4 * mp.pi * elliptic_square(angle) - mp.pi ** 3) / mp.sin(angle) ** 2
subtracted_integral = mp.quad(
subtracted, [0, mp.pi / 8, mp.pi / 4], method='gauss-legendre')
gamma_constant = mp.pi ** 3 - (
mp.gamma(mp.mpf(1) / 4) ** 4 + 16 * mp.gamma(mp.mpf(3) / 4) ** 4) / 8
fourth_second = (-2 * fourth_zero + 3 * subtracted_integral
+ 3 * mp.pi ** 2 - 3 * gamma_constant) / 16
def complementary_product(variable):
if not variable:
return mp.mpf(0)
angle = variable ** 2
sine_squared = mp.sin(angle) ** 2
cosine_squared = mp.cos(angle) ** 2
return (4 * variable * mp.elliprf(0, cosine_squared, 1)
* mp.elliprf(0, sine_squared, 1))
mixed_zero = mp.quad(complementary_product,
[0, mp.sqrt(mp.pi / 8), mp.sqrt(mp.pi / 4)],
method='tanh-sinh')
result = dict(zip(EVEN_KEYS, (third_zero, third_second, fourth_zero,
fourth_second, mixed_zero)))
if not all(mp.isfinite(value) and value > 0 for value in result.values()):
raise ArithmeticError('Nonfinite or nonpositive independent value')
return result
def check(rows):
entries = {tuple(int(row['params'][name]) for name in ('a', 'b', 'k')): row
for row in rows}
if not set(EVEN_KEYS).issubset(entries):
raise ValueError('Missing even moments')
lower = reference_values(140)
higher = reference_values(180)
records = []
with mp.workdps(210):
for key in EVEN_KEYS:
reference = higher[key]
agreement = abs(lower[key] - reference) / abs(reference)
produced = mp.mpf(entries[key]['at_200'])
difference = abs(produced - reference) / abs(reference)
written = entries[key]['number']
text_center = mp.mpf(written)
decimal_exponent = Decimal(written).as_tuple().exponent
written_error_ulps = abs(text_center - reference) / mp.mpf(10) ** decimal_exponent
if not agreement < mp.mpf('1e-125') or not difference < mp.mpf('1e-125'):
raise ArithmeticError(('Independent check failed', key, agreement, difference))
if not written_error_ulps < mp.mpf('0.6'):
raise ArithmeticError(('Stored last-place claim failed', key, written_error_ulps))
records.append({
'params': key, 'independent_at_140': mp.nstr(lower[key], 140),
'independent_at_180': mp.nstr(reference, 180),
'independent_precision_agreement': mp.nstr(agreement, 12),
'relative_difference_from_bessel_quadrature': mp.nstr(difference, 12),
'stored_error_ulps': mp.nstr(written_error_ulps, 12),
'finite_and_positive': True})
return records
if __name__ == '__main__':
rows = [json.loads(line) for line in Path(sys.argv[1]).read_text().splitlines() if line]
print(json.dumps({'independent_checks': check(rows), 'complete': True}, indent=2))