import numberdb.sage as numberdb
from sage.rings.rational_field import QQ
from sage.symbolic.ring import SR
import sage.rings.polynomial.laurent_polynomial_ring
from sage.combinat.sf.sf import SymmetricFunctions
from sage.combinat.partition import Partitions
from math import factorial, prod
from collections import Counter
def orthogonal_weingarten(N, mu):
k = sum(mu)
sf = SymmetricFunctions(QQ)
p, s, J = sf.power(), sf.schur(), sf.jack(t=QQ(2)).J()
z = prod(r**m * factorial(m) for r, m in Counter(mu).items())
total = QQ(0)
for lam in Partitions(k):
if len(lam) > N:
continue
f = s[tuple(2*r for r in lam)].scalar(p[[1]*(2*k)])
c = p(J[lam]).coefficient(mu)
d = prod(N + 2*j - i for i, row in enumerate(lam) for j in range(row))
total += QQ(f*c) / QQ(d)
return QQ(2**len(mu)*z) * total / QQ(factorial(2*k))
print(orthogonal_weingarten(13, (3, 1)))
The attached generator evaluates the spherical expansion over $\mathbb Q$.
It obtains zonal spherical coefficients from the expansion in power-sum symmetric polynomials of Jack polynomials $J_\lambda^{(2)}$, with the normalization of [2], equations (3.5)-(3.7), and computes $f^{2\lambda}$ by the hook-length formula. The values are exact fractions, without numerical rounding. All 242 coefficients were independently recovered from pairing Gram inverse equations by rational Gaussian elimination, using up to 10,395 pairings. The checks also verify 57 moments $\mathbb E R_{11}^{2k}=(2k-1)!!/\prod_{j=0}^{k-1}(N+2j)$ from the uniform distribution on the unit sphere, and 64 explicit coefficients through $k=3$ from [2], section 7. Nonzero controls and an altered coefficient test the checking procedure. The reader program agrees with every entry; the zonal specializations and the singular $N=1$ convention are checked separately in the attachments.