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 collections import Counter
from math import factorial, prod


def symplectic_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
        doubled_rows = tuple(r for r in lam for _ in (0, 1))
        f = s[doubled_rows].scalar(p[[1]*(2*k)])
        c = p(J[lam.conjugate()]).coefficient(mu)
        d = prod(2*N - 2*i + j 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(symplectic_weingarten(2, (2,)))  # 1/40