generate.py

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

3325 bytes, as of the version from 2026-10-26 23:52 (current). Recorded here, not run.

"""Airy derivative zeta values -- numberdb.org/T532.

Run with SageMath or a modular Sage installation:
    sage -pip install numberdb
    sage -python generate.py             # verify every stored entry
    sage -python checks.py               # independent numerical checks
    sage -python generate.py --fill      # write entries to this draft

Writing requires NUMBERDB_API_KEY, or NUMBERDB_KEY_FROM_STDIN=1 with the key
on stdin. --fill invokes the client's entry writer; it does not make a draft
public. The table is left for a separate critique and human publication.

Let t = 1/kappa = Gamma(1/3)/(3^(1/3) Gamma(2/3)). The exact polynomials
P_2(t)=t, P_3(t)=1 and
P_n(t)=sum(P_j(t)*P_(n-j)(t), j=2..n-2)/(n-2)
follow from h'=1-z*h² for h=Ai/Ai'. The canonical product gives
Z_N(n)=(-1)^(n-1) [z^(n-2)] h(z) = P_n(t).
All polynomial coefficients and divisions are rational. Only the final
substitution uses Arb balls. A 96-bit guard is tested in checks.py; the
100-digit default uses 429 bits. Exactly constant polynomials stay exact.

Sources: Voros, math-ph/9811001, Section 3, especially (28),(29),(43.2-3);
NIST DLMF 9.2.1, 9.2.3, 9.2.4 and 9.9. See argument-selection.md.
"""
import os
import re
import sys
from functools import lru_cache
import numberdb.sage as numberdb
from sage.rings.rational_field import QQ
from sage.rings.polynomial.polynomial_ring_constructor import PolynomialRing
from sage.rings.real_arb import RealBallField

MAX_ORDER = 20
WORKING_GUARD = 96
R = PolynomialRing(QQ, 't')

@lru_cache(maxsize=None)
def polynomials(max_order):
    p = [R(0), R(0), R.gen(), R(1)]
    for n in range(4, max_order + 1):
        p.append(sum((p[j] * p[n-j] for j in range(2, n-1)), R(0)) / QQ(n-2))
    return tuple(p)

@lru_cache(maxsize=None)
def inverse_kappa(bits):
    b = RealBallField(bits)
    return b(QQ(1)/3).gamma() / (b(3)**(QQ(1)/3) * b(QQ(2)/3).gamma())

class AiryDerivativeZeta(numberdb.Generator):
    table = 'T532'
    parameters = ('n',)
    type = 'R'
    digits = 100
    rigour = 'proven'

    def enumerate(self, max_order=MAX_ORDER):
        for n in range(2, max_order + 1):
            yield {'n': str(n)}

    def value(self, params, digits):
        n = int(params['n'])
        if n < 2:
            raise ValueError('The convergent family has integer n >= 2')
        poly = polynomials(max(n, MAX_ORDER))[n]
        if poly.degree() == 0:
            return {'number': QQ(poly[0]), 'comment': r'The inverse cubes sum to $1$. CITE{Voros}'}
        value = poly(inverse_kappa(numberdb.bits(digits, losing=WORKING_GUARD)))
        assert value.is_finite() and value > 0
        return value

def read_key_if_requested():
    if os.environ.get('NUMBERDB_KEY_FROM_STDIN') == '1':
        key = sys.stdin.read().strip()
        if re.match(r'^[A-Z_][A-Z0-9_]*\s*=', key):
            key = key.split('=', 1)[1].strip().strip('\"\'')
        os.environ['NUMBERDB_API_KEY'] = key

if __name__ == '__main__':
    read_key_if_requested()
    generator = AiryDerivativeZeta()
    if '--fill' in sys.argv or os.environ.get('NUMBERDB_PUBLISH') == '1':
        print(generator.publish(message='Compute Neumann Airy moments with exact polynomials and Arb enclosures'))
    else:
        report = generator.verify(sample=None)
        print(report)
        sys.exit(0 if report.ok else 1)