generate.py

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

5070 bytes, as of the version from 2026-11-01 14:17. Recorded here, not run.

"""Fourier coefficients of monstrous McKay–Thompson series -- numberdb.org/T652.

Run with SageMath and the attached seeds.json and class-data.json:
    sage -pip install numberdb
    sage -python generate.py              # verify all stored entries
    sage -python generate.py --fill       # write entries to the existing draft
The --fill mode needs NUMBERDB_API_KEY. It never makes the draft public.

Borcherds, Monstrous moonshine and monstrous Lie superalgebras (1992),
§9, equation (9.1), determines the coefficients from n=1,2,3,5 and the
squared-class series. The seed values were transcribed from Conway–Norton,
Table 4, printed p. 334; H_0 in that source is not a Fourier constant term.
All divisions by two below require exact divisibility. Class power maps are
frozen from GAP CTblLib 1.3.11. Identity coefficients are computed internally
from E4^3/Delta-744, but excluded from enumeration because T79 covers them.
See argument-selection.md for both discrete indices and the finite bounds.
"""
import numberdb.sage as numberdb
from sage.rings.integer_ring import ZZ
from functools import lru_cache
from pathlib import Path
import json
import os
import sys

HERE = Path(__file__).resolve().parent
SEEDS = json.loads((HERE / 'seeds.json').read_text())
CLASS_DATA = json.loads((HERE / 'class-data.json').read_text())
CLASSES = tuple(g for g, order in zip(CLASS_DATA['class_names'],
                                    CLASS_DATA['class_orders'])
                if 1 < order <= 10)
SQUARE = CLASS_DATA['power_maps']['2']
assert len(CLASSES) == 32
assert set(SQUARE) == set(CLASSES) | {'1A'}
assert set(SQUARE.values()) <= set(SQUARE)


def product(a, b, degree):
    return [sum((a[j] * b[n-j] for j in range(max(0,n-len(b)+1),
                min(n,len(a)-1)+1)), ZZ(0)) for n in range(degree+1)]


def identity_series(nmax):
    """Exact E4^3 / product(1-q^m)^24, shifted by q^-1 and minus 744."""
    degree = nmax + 1
    e4 = [ZZ(1)] + [ZZ(240)*sum(ZZ(d)**3 for d in range(1,n+1)
                               if n % d == 0) for n in range(1,degree+1)]
    numerator = product(product(e4,e4,degree),e4,degree)
    denominator = [ZZ(1)] + [ZZ(0)]*degree
    for m in range(1,degree+1):
        for _ in range(24):
            for k in range(degree,m-1,-1):
                denominator[k] -= denominator[k-m]
    quotient = []
    for n in range(degree+1):
        quotient.append(numerator[n] - sum(denominator[j]*quotient[n-j]
                                           for j in range(1,n+1)))
    assert quotient[0] == 1 and quotient[1] == 744
    return {n:quotient[n+1] for n in range(1,nmax+1)}


def half(value):
    value = ZZ(value)
    assert value % 2 == 0
    return value // 2


@lru_cache(maxsize=None)
def coefficients(nmax=30):
    assert nmax >= 5
    data = {g:{int(n):ZZ(v) for n,v in SEEDS[g].items()} for g in CLASSES}
    data['1A'] = identity_series(nmax)
    for n in range(4,nmax+1):
        if n == 5:
            continue
        k, r = divmod(n,4)
        for g in CLASSES:
            a = data[g]
            b = data[SQUARE[g]]
            if r == 0:
                value = a[2*k+1] + half(a[k]**2-b[k])
                value += sum(a[j]*a[2*k-j] for j in range(1,k))
            elif r == 1:
                value = a[2*k+3]-a[2]*a[2*k]+half(a[2*k]**2+b[2*k])
                value += half(a[k+1]**2-b[k+1])
                value += sum(a[j]*a[2*k-j+2] for j in range(1,k+1))
                value += sum(b[j]*a[4*k-4*j] for j in range(1,k))
                value += sum((-1)**j*a[j]*a[4*k-j] for j in range(1,2*k))
            elif r == 2:
                value = a[2*k+2]+sum(a[j]*a[2*k-j+1] for j in range(1,k+1))
            else:
                value = a[2*k+4]-a[2]*a[2*k+1]-half(a[2*k+1]**2-b[2*k+1])
                value += sum(a[j]*a[2*k-j+3] for j in range(1,k+2))
                value += sum(b[j]*a[4*k-4*j+2] for j in range(1,k+1))
                value += sum((-1)**j*a[j]*a[4*k-j+2] for j in range(1,2*k+1))
            a[n] = ZZ(value)
    return data


class McKayThompsonCoefficients(numberdb.Generator):
    table = 'T652'
    parameters = ('g','n')
    type = 'Z'
    rigour = 'exact'
    files = ('generate.py','seeds.json','class-data.json','argument-selection.md',
             'checks.py','check.py','published-controls.json','gap_export.py')

    def enumerate(self, nmax=30):
        for g in CLASSES:
            for n in range(1,nmax+1):
                yield {'g':g, 'n':n}

    def value(self, params, digits):
        n = int(params['n'])
        return coefficients(max(30,n))[params['g']][n]


if __name__ == '__main__':
    if os.environ.get('NUMBERDB_KEY_FROM_STDIN') == '1':
        os.environ['NUMBERDB_API_KEY'] = sys.stdin.read().strip()
    generator = McKayThompsonCoefficients()
    if '--fill' in sys.argv or os.environ.get('NUMBERDB_PUBLISH') == '1':
        print(generator.publish(message='Fill 960 signed exact coefficients after independent checks'))
    else:
        report = generator.verify(sample=None)
        print(report)
        sys.exit(0 if report.ok else 1)