back to table · edit · history · where entries came from · files · download
4046 bytes, as of the version from 2026-11-09 09:52 (current). Recorded here, not run.
"""Multiple Mahler measures with a shared phase -- numberdb.org/T685.
Run with SageMath and the numberdb package:
sage -pip install numberdb
sage -python generate.py # verify the stored table
sage -python generate.py --submit # fill your draft using your API key
This computes n=3,4 and every 1 <= k <= 12. Parameters are (n,k), where n
counts all summands and k counts logarithmic factors. One phase is shared;
all other phases are private. See argument-selection.md and certification.md.
Three-term values use certified Arb integration, with analytic log branches.
Four-term values integrate a polynomial truncation of the conditional mean,
with a proved uniform remainder bound. WORKING_GUARD covers arithmetic and
integration losses; every output must retain at least 24 extra accurate bits.
No stored digits or data files enter the computation.
"""
import os
import sys
from functools import lru_cache
import numberdb.sage as numberdb
from sage.rings.integer_ring import ZZ
from sage.rings.rational_field import QQ
from sage.rings.real_arb import RealBallField
from sage.rings.complex_arb import ComplexBallField
from sage.rings.polynomial.polynomial_ring_constructor import PolynomialRing
WORKING_GUARD = 112
def field_bits(digits):
return numberdb.bits(digits, losing=WORKING_GUARD)
@lru_cache(maxsize=None)
def four_term_polynomial(digits):
bits = field_bits(digits)
R = RealBallField(bits)
P = PolynomialRing(R, 'u')
# 4^(-N) is smaller than working rounding; the bound includes all tails.
N = (bits + 25) // 2
coefficients = [R(1)] + [
-2 * R(2*m).zeta() / ((2*m+1) * ZZ(4)**m)
for m in range(1, N+1)
]
h = P(coefficients)
error = R(QQ(16) / (3*(2*N+3)*ZZ(4)**(N+1)))
return h, error, N
@lru_cache(maxsize=None)
def moment(n, k, digits=100):
n, k = ZZ(n), ZZ(k)
if n not in (3, 4) or k < 1:
raise ValueError('Implemented summand counts are n=3,4; require k>=1.')
bits = field_bits(digits)
R = RealBallField(bits)
if n == 3:
C = ComplexBallField(bits)
pi = C.pi()
integral = C.integral(
lambda t, analytic: (2*(pi*t/3).cos()).log(analytic=analytic)**k,
0, 1, rel_tol=R(2)**(-bits+24), abs_tol=R(2)**(-bits+24))
if not integral.imag().contains_zero():
raise ArithmeticError('Real integral has nonzero imaginary enclosure.')
result = (2*integral/3).real()
else:
h, error, _ = four_term_polynomial(digits)
power = h**k
result = sum((c / (k+2*j+1) for j,c in enumerate(power.list())), R(0))
result = result.add_error(k*(1+error)**(k-1)*error)
if not result.is_finite() or result.accuracy() < numberdb.bits(digits)+24:
raise ArithmeticError('Insufficient finite relative accuracy: %s' % result)
return result
class MultipleMahler(numberdb.Generator):
table = 'T685'
parameters = ('n', 'k')
type = 'R'
digits = 100
rigour = 'proven'
files = ('generate.py', 'argument-selection.md', 'certification.md',
'checks.py', 'numerical-evidence.json', 'skill-manifest.json')
def enumerate(self):
for n in (3, 4):
for k in range(1, 13):
yield {'n': str(n), 'k': str(k)}
def value(self, params, digits):
return moment(int(params['n']), int(params['k']), digits)
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].strip().isidentifier():
key = key.split('=',1)[1].strip().strip('\"\'')
os.environ['NUMBERDB_API_KEY'] = key
generator = MultipleMahler()
if '--submit' in sys.argv or os.environ.get('NUMBERDB_PUBLISH') == '1':
print(generator.publish(message='Fill the shared-phase multiple Mahler family with certified values'))
else:
report = generator.verify(sample=None)
print(report)
sys.exit(0 if report.ok else 1)