"""Sharp lower constants in the Khintchine inequality -- numberdb.org/T605.

Reproduce with SageMath (or passagemath with numberdb installed):
    sage -python generate.py                 # verify every stored entry
    sage -python generate.py --fill          # write entries to this draft
    sage -python checks.py                   # independent numerical checks
Set NUMBERDB_API_KEY for authenticated draft access and writes. The --fill
operation writes numbers; it cannot offer, review, or publish the draft.

The definition uses real coefficients, independent equiprobable signs, and
probability expectation. A_p is the largest universal lower norm constant.
Mordhorst (2017), Theorem 1.1, gives min(1, 2^(1/2-1/p), G(p)).
All powers are positive real. Exact rationals are formed before ball conversion.
Strict Arb comparisons select the active branch. An unresolved comparison
raises; it never chooses by midpoints. p >= 2 is exactly 1 by norm monotonicity.
Integer exponents in the power-of-two branch are evaluated in QQ.

The 96-bit guard is measured in numerical-evidence.json. The initial selection
and every cutoff are explained in argument-selection.md; no external input
file or source code controls enumeration. See skill-manifest.json for the
archived local/live instruction hashes applied to this run.
"""
import os
import sys
import numberdb.sage as numberdb
from sage.rings.rational_field import QQ
from sage.rings.real_arb import RealBallField

WORKING_GUARD = 96
MAX_REFERENCE_DENOMINATOR = 6
MAX_APPLICATION_INDEX = 20


def selection_groups():
    intrinsic = {QQ(1), QQ(2)}
    application = {QQ(2*m) / QQ(m+1) for m in range(2, MAX_APPLICATION_INDEX+1)}
    reference = {QQ(a)/QQ(b)
                 for b in range(1, MAX_REFERENCE_DENOMINATOR+1)
                 for a in range(1, 2*b)} - intrinsic - application
    return {'intrinsic': sorted(intrinsic),
            'documented_application': sorted(application),
            'reference_sample': sorted(reference)}


def evaluate(p, digits=100):
    p = QQ(p)
    if p <= 0:
        raise ValueError('p must be positive')
    if p >= 2:
        return QQ(1), 'unit'
    field = RealBallField(numberdb.bits(digits, losing=WORKING_GUARD))
    exponent = QQ(1)/QQ(2) - QQ(1)/p
    power = (field(2).log() * field(exponent)).exp()
    gamma = (field(2).log()/2 +
             (field((p+1)/2).gamma().log() - field.pi().log()/2)/field(p)).exp()
    assert power.is_finite() and gamma.is_finite()
    if power < gamma and power < 1:
        if exponent.denominator() == 1:
            return QQ(2)**int(exponent), 'two-point-exact'
        return power, 'two-point'
    if gamma < power and gamma < 1:
        return gamma, 'gaussian'
    raise ArithmeticError('Ball comparison did not separate the branches at p=%s' % p)


class KhintchineLower(numberdb.Generator):
    table = 'T605'
    parameters = ('p',)
    type = 'R'
    digits = 100
    rigour = 'proven'
    files = ('generate.py', 'checks.py', 'argument-selection.md',
             'numerical-evidence.json', 'skill-manifest.json', 'sources.md')

    def enumerate(self):
        groups = selection_groups()
        for p in sorted(set().union(*map(set, groups.values()))):
            yield {'p': str(p)}

    def value(self, params, digits):
        p = QQ(params['p'])
        value, branch = evaluate(p, digits)
        comment = None
        if p == 1:
            comment = r'$A_1=1/\sqrt{2}$, the sharp first-moment constant.'
        elif p == 2:
            comment = r'$A_2=1$ follows from $\mathbb E(\sum_j a_j\varepsilon_j)^2=\sum_j a_j^2$.'
        elif branch == 'two-point-exact':
            comment = r'The exponent $1/2-1/p$ is an integer, so this value is exact.'
        return {'number': value, 'comment': comment} if comment else value


if __name__ == '__main__':
    if os.environ.get('NUMBERDB_KEY_FROM_STDIN') == '1':
        os.environ['NUMBERDB_API_KEY'] = sys.stdin.read().strip()
    generator = KhintchineLower()
    if '--fill' in sys.argv or os.environ.get('NUMBERDB_PUBLISH') == '1':
        print(generator.publish(message='Fill independently checked sharp lower Khintchine constants'))
    else:
        report = generator.verify(sample=None)
        print(report)
        sys.exit(0 if report.ok else 1)
