back to table · edit · history · where entries came from · files · download
10538 bytes, as of the version from 2026-10-27 01:03 (current). Recorded here, not run.
"""Independent numerical and identity checks, with decimal-interval semantics.
python -m pip install numberdb python-flint==0.6.0 mpmath==1.3.0
python checks.py T351 --live
No API writes. --document checks a local whole-document JSON file instead.
"""
import argparse
import json
import math
from decimal import Decimal
from fractions import Fraction
from pathlib import Path
import mpmath as mp
import numberdb
from flint import arb
from values import PARAMETERS, bernoulli, endpoint, evaluate, precision, rational, selection
def flatten(document):
result = {}
depth = len(document["Parameters"])
def descend(node, keys):
if len(keys) == depth:
if isinstance(node, dict):
if "number" not in node:
return
node = node["number"]
result[tuple(keys)] = str(node)
elif isinstance(node, dict):
for key, value in node.items():
if key not in ("equals", "comment", "param-latex", "proof", "number"):
descend(value, keys + [key])
descend(document["Numbers"], [])
return result
def bounds(text):
midpoint = Fraction(text)
if "." not in text and "e" not in text.lower():
return midpoint, midpoint, Fraction(0)
unit = Fraction(10) ** Decimal(text).as_tuple().exponent
return midpoint - unit, midpoint + unit, unit
def numeric_fraction(value):
return mp.mpf(value.numerator) / value.denominator
def mp_argument(identity):
if identity == "sqrt(2)":
return mp.sqrt(2)
if identity == "1/sqrt(2)":
return 1 / mp.sqrt(2)
return numeric_fraction(Fraction(identity))
def mp_value(tid, params):
arguments = {key: mp_argument(value) for key, value in params.items()}
if tid == "T14":
return mp.zeta(arguments["s"])
if tid == "T94":
return mp.zeta(arguments["s"], arguments["x"])
if tid in ("T228", "T246"):
value = mp.diff(mp.zeta, arguments["s"])
return value if tid == "T228" else value / mp.zeta(arguments["s"])
if tid == "T51":
return mp.agm(arguments["a"], arguments["b"])
if tid == "T176":
return mp.polylog(arguments["s"], arguments["x"])
if tid == "T178":
return mp.hyp2f1(1 / arguments["r"], 1 - 1 / arguments["r"], 1, arguments["x"])
if tid == "T351":
return mp.gammainc(arguments["a"], arguments["x"], mp.inf)
if tid == "T352":
return mp.gammainc(arguments["a"], 0, arguments["x"])
if tid == "T353":
return mp.betainc(arguments["a"], arguments["b"], 0, arguments["x"])
raise ValueError("unknown table")
def upper_half_integer(shape, point):
count = int(2 * shape)
if count % 2 == 0:
integer = count // 2
return math.factorial(integer - 1) * (-point).exp() * sum(
point ** power / math.factorial(power) for power in range(integer))
value = arb.pi().sqrt() * point.sqrt().erfc()
current = Fraction(1, 2)
while current < shape:
value = rational(current) * value + point ** rational(current) * (-point).exp()
current += 1
return value
def identity_value(tid, params):
if tid in ("T351", "T352"):
shape = Fraction(params["a"])
upper = upper_half_integer(shape, rational(params["x"]))
return upper if tid == "T351" else rational(shape).gamma() - upper
if tid == "T353":
shape = rational(params["a"])
return shape.gamma() ** 2 / (2 * (2 * shape).gamma())
if tid == "T178":
signature = int(params["r"])
first = rational(Fraction(signature + 1, 2 * signature))
second = rational(Fraction(2 * signature - 1, 2 * signature))
return arb.pi().sqrt() / (first.gamma() * second.gamma())
if tid == "T51":
base = 2 * arb.pi() ** rational(Fraction(3, 2)) / rational("1/4").gamma() ** 2
return base * arb(2).sqrt() if params["b"] == "sqrt(2)" else base
if tid == "T176":
order = int(params["s"])
if params["x"] == "-1":
return -(1 - arb(2) ** (1 - order)) * arb(order).zeta()
logarithm = arb(2).log()
if params["x"] == "1/2" and order == 2:
return arb.pi() ** 2 / 12 - logarithm ** 2 / 2
if params["x"] == "1/2" and order == 3:
return 7 * arb(3).zeta() / 8 - arb.pi() ** 2 * logarithm / 12 + logarithm ** 3 / 6
if params["x"] == "1/3" and order == 2:
return -rational("-1/2").polylog(2) - rational("2/3").log() ** 2 / 2
if params["x"] == "1/9" and order == 2:
return 6 * rational("1/3").polylog(2) - arb.pi() ** 2 / 3 + arb(3).log() ** 2
if tid == "T14":
index = Fraction(params["s"])
if index.denominator == 1 and index >= 2 and index.numerator % 2 == 0:
integer = index.numerator
return ((-1) ** (integer // 2 + 1) * rational(bernoulli(integer)) *
(2 * arb.pi()) ** integer / (2 * math.factorial(integer)))
if tid == "T94" and params["x"] == "1/2":
index = rational(params["s"])
return (arb(2) ** index - 1) * index.zeta()
if tid == "T228":
index = Fraction(params["s"])
if index == 0:
return -(2 * arb.pi()).log() / 2
if index.denominator == 1 and index < 0 and index.numerator % 2 == 0:
degree = -index.numerator
return ((-1) ** (degree // 2) * math.factorial(degree) * arb(degree + 1).zeta()
/ (2 * (2 * arb.pi()) ** degree))
if tid == "T246":
if params["s"] == "0":
return (2 * arb.pi()).log()
if params["s"] == "1/2":
return (arb.pi().log() - rational("1/4").digamma()) / 2
return None
def check(tid, document):
actual = flatten(document)
expected = {tuple(params[name] for name in PARAMETERS[tid]) for params in selection(tid)}
if set(actual) != expected:
raise AssertionError("enumeration differs from stored rows")
records = []
with precision(100):
control = upper_half_integer(Fraction(1), arb(1))
if not control.is_finite() or not (control - arb(-1).exp()).contains(0):
raise AssertionError("nonzero control failed")
for params in selection(tid):
key = tuple(params[name] for name in PARAMETERS[tid])
text = actual[key]
lower, upper, unit = bounds(text)
row = {"parameters": params, "number": text}
with precision(100):
value = evaluate(tid, params)
if isinstance(value, Fraction):
if lower != upper or value != lower:
raise AssertionError("exact arithmetic mismatch: " + str(params))
row["exact"] = True
else:
if not value.is_finite() or value.rel_accuracy_bits() < 360:
raise AssertionError("nonfinite or uninformative result")
measured_lower, measured_upper = endpoint(value.lower()), endpoint(value.upper())
if not lower < measured_lower <= measured_upper < upper:
raise AssertionError("stored decimal fails containment: " + str(params))
if lower + 10 * unit < measured_lower <= measured_upper < upper + 10 * unit:
raise AssertionError("negative control accepted a shifted value")
row["certified_accuracy_bits"] = value.rel_accuracy_bits()
row["decimal_contains_enclosure"] = True
row["shifted_decimal_rejected"] = True
identity = identity_value(tid, params)
if identity is not None:
measured = rational(value) if isinstance(value, Fraction) else value
if not identity.is_finite():
raise AssertionError("nonfinite identity check")
if identity != 0 and identity.rel_accuracy_bits() < 360:
raise AssertionError("wide identity check")
if not (measured - identity).contains(0):
raise AssertionError("identity disagreement")
row["independent_identity"] = True
if unit and not lower < endpoint(identity.lower()) <= endpoint(identity.upper()) < upper:
raise AssertionError("independent identity enclosure outside decimal")
previous = None
for working in (150, 200):
with mp.workdps(working):
independent = mp_value(tid, params)
if not mp.isfinite(independent):
raise AssertionError("nonfinite mpmath value")
if unit:
if not numeric_fraction(lower) < independent < numeric_fraction(upper):
raise AssertionError("mpmath disagreement: " + str(params))
elif abs(independent - numeric_fraction(lower)) > mp.mpf(10) ** -140 * max(1, abs(independent)):
raise AssertionError("mpmath exact-value disagreement")
if previous is not None and abs(independent - mp.mpf(previous)) > mp.mpf(10) ** -139 * max(1, abs(independent)):
raise AssertionError("independent precisions disagree")
previous = mp.nstr(independent, 145)
row["mpmath_150_200_pass"] = True
records.append(row)
return {"tid": tid, "entries": len(records), "python_flint": "0.6.0",
"mpmath": "1.3.0", "working_bits": numberdb.bits(100, losing=160),
"minimum_certified_accuracy_bits": min(
row.get("certified_accuracy_bits", 10 ** 6) for row in records),
"identity_checks": sum(row.get("independent_identity", False) for row in records),
"exact_rows": sum(row.get("exact", False) for row in records),
"nonzero_control_passed": True, "records": records}
if __name__ == "__main__":
parser = argparse.ArgumentParser(description=__doc__)
parser.add_argument("tid", choices=PARAMETERS)
parser.add_argument("--document", type=Path)
parser.add_argument("--live", action="store_true")
parser.add_argument("--output", type=Path)
arguments = parser.parse_args()
document = (json.loads(arguments.document.read_text()) if arguments.document
else numberdb.table(arguments.tid))
evidence = check(arguments.tid, document)
if arguments.output:
arguments.output.write_text(json.dumps(evidence, indent=2) + "\n")
print(json.dumps({key: value for key, value in evidence.items() if key != "records"}))