back to table · edit · history · where entries came from · files · download
5545 bytes, as of the version from 2026-10-28 08:03 (current). Recorded here, not run.
"""Jack polynomials -- numberdb.org/T575
Run with SageMath:
sage -pip install numberdb
sage -python generate.py # verify the stored draft; needs its key
sage -python generate.py --fill # fill entries; never publishes the table
The table stores the integral J form over ZZ[x1,...,xn,a], with a = alpha.
The two entry parameters are n and lambda, in that order. The default range
is documented in argument-selection.md. All arithmetic is exact.
Construction: express monomial symmetric functions in power sums by the
Moebius function of set partitions, then perform Gram-Schmidt over QQ(a)
with <p_mu,p_nu> = delta_mu,nu z_mu a^length(mu). Multiply the monic result
by product_s (a arm(s)+leg(s)+1). Only then specialize to n alphabet
variables. This avoids depending on Sage's symmetric-function implementation.
"""
import numberdb.sage as numberdb
from sage.rings.rational_field import QQ
from sage.rings.integer_ring import ZZ
from sage.rings.polynomial.polynomial_ring_constructor import PolynomialRing
from collections import Counter
from functools import lru_cache
from math import factorial
import os
import sys
A = PolynomialRing(QQ, 'a')
K = A.fraction_field()
a = K(A.gen())
# Complete degree slices, chosen from measured expansion lengths.
DEGREE_BOUNDS = ((2, 5), (3, 5), (4, 4))
def partitions(d, maximum=None):
"""Partitions in decreasing lexicographic order, including the empty one."""
if d == 0:
yield ()
return
for first in range(min(d, d if maximum is None else maximum), 0, -1):
for tail in partitions(d-first, first):
yield (first,) + tail
def set_partitions(items):
if not items:
yield ()
return
first, *tail = items
for blocks in set_partitions(tail):
yield ((first,),) + blocks
for i in range(len(blocks)):
yield blocks[:i] + ((first,) + blocks[i],) + blocks[i+1:]
def z_value(shape):
z = 1
for part, multiplicity in Counter(shape).items():
z *= part**multiplicity * factorial(multiplicity)
return ZZ(z)
def monomial_in_p(shape):
divisor = 1
for count in Counter(shape).values():
divisor *= factorial(count)
out = {}
for blocks in set_partitions(list(range(len(shape)))):
mu = tuple(sorted((sum(shape[i] for i in block) for block in blocks), reverse=True))
c = QQ(1) / QQ(divisor)
for block in blocks:
c *= (-1)**(len(block)-1) * factorial(len(block)-1)
out[mu] = out.get(mu, K(0)) + c
return {mu: K(c) for mu, c in out.items() if c}
def pairing(left, right):
return sum((c * right.get(mu, K(0)) * z_value(mu) * a**len(mu)
for mu, c in left.items()), K(0))
def hook_product(shape):
out = K(1)
for i, width in enumerate(shape):
for j in range(width):
arm = width-j-1
leg = sum(other > j for other in shape[i+1:])
out *= a*arm+leg+1
return out
@lru_cache(maxsize=None)
def integral_basis(d):
orthogonal = {}
norms = {}
for shape in reversed(list(partitions(d))):
v = monomial_in_p(shape)
for mu, q in orthogonal.items():
factor = pairing(v, q) / norms[mu]
for rho, c in q.items():
v[rho] = v.get(rho, K(0)) - factor*c
v = {rho: c for rho, c in v.items() if c}
orthogonal[shape] = v
norms[shape] = pairing(v, v)
assert norms[shape] != 0
return {shape: {mu: c*hook_product(shape) for mu, c in v.items()}
for shape, v in orthogonal.items()}
def jack(n, shape):
if n < 1 or len(shape) > n or any(k < 1 for k in shape) or tuple(sorted(shape, reverse=True)) != shape:
raise ValueError('n must be positive and lambda a partition of length at most n')
names = ['x%d' % i for i in range(1, n+1)]
X = PolynomialRing(K, names)
xs = X.gens()
power = {k: sum((x**k for x in xs), X(0)) for k in range(1, sum(shape)+1)}
f = X(0)
for mu, coeff in integral_basis(sum(shape))[shape].items():
term = X(coeff)
for k in mu:
term *= power[k]
f += term
R = PolynomialRing(ZZ, names+['a'])
terms = {}
for exponents, coeff in f.dict().items():
if n == 1:
exponents = (int(exponents),)
# Coercion refuses any non-polynomial denominator, then any noninteger.
for power_a, c in A(coeff).dict().items():
terms[tuple(exponents)+(int(power_a),)] = ZZ(c)
return R(terms)
class JackPolynomials(numberdb.Generator):
table = 'T575'
parameters = ('n', 'lambda')
type = 'Z[]'
rigour = 'exact'
def enumerate(self):
for n, upper_degree in DEGREE_BOUNDS:
for d in range(2, upper_degree+1):
for shape in partitions(d):
if len(shape) <= n:
yield {'n': str(n), 'lambda': ','.join(map(str, shape))}
def value(self, params, digits):
return jack(int(params['n']), tuple(int(s) for s in params['lambda'].split(',')))
if __name__ == '__main__':
if os.environ.get('NUMBERDB_KEY_FROM_STDIN'):
numberdb.configure(api_key=sys.stdin.read().strip())
generator = JackPolynomials()
if '--fill' in sys.argv or os.environ.get('NUMBERDB_PUBLISH'):
print(generator.publish(message='Fill exact integral Jack polynomials after independent checks'))
else:
report = generator.verify(sample=None)
print(report)
sys.exit(0 if report.ok else 1)