OFFSET
0,3
LINKS
Philip Turecek, Table of n, a(n) for n = 0..8
FORMULA
a(n) = Sum_{1*s_1+2*s_2+...=n, 1*t_1+2*t_2+...=n, 1*u_1+2*u_2+...=n, 1*v_1+2*v_2+...=n} (fixA[s_1, s_2, ...;t_1, t_2, ...;u_1, u_2, ...;v_1, v_2, ...]/ (1^s_1*s_1!*2^s_2*s_2!*... *1^t_1*t_1!*2^t_2*t_2!*... *1^u_1*u_1!*2^u_2*u_2!*... *1^v_1*v_1!*2^v_2*v_2!*...)) where fixA[...] = Product_{i, j, k>=1} ( (Sum_{d|lcm(i, j, k)} (d*v_d))^(s_i*t_j*u_k *lcm(i, j, k)/(i*j*k))).
a(n) asymptotic to n^(n^3)/(n!^4).
PROG
(Sage)
Pol.<x> = InfinitePolynomialRing(QQ)
@cached_function
def Z(n):
if n == 0: return Pol.one()
return sum(x[k]*Z(n-k) for k in (1..n))/n
@cached_function
def monprod(M):
p = Pol.one()
V = [m.variables() for m in M]
T = cartesian_product(V)
for t in T:
r = [Pol.varname_key(str(u))[1] for u in t]
j = [Pol(M[u[0]]).degree(u[1]) for u in enumerate(t)]
lcm_r = lcm(r)
p *= x[lcm_r]^(prod(r)/lcm_r*prod(j))
return p
@cached_function
def pol_isotop(n, k):
P = Z(n)
p = Pol.zero()
coeffs = P.coefficients()
mons = P.monomials()
C = cartesian_product(k*[mons])
Csorted = [tuple(sorted(u)) for u in C]
Cset = set(Csorted)
for c in Cset:
p += Csorted.count(c)*prod([coeffs[mons.index(u)] for u in c])*monprod(c)
return p
@cached_function
def rule_sub(r, m):
D = 0
for d in divisors(r):
try: D += d*m.degrees()[-d-1]
except: break
return D
def a(n, k=3):
P = Z(n)
coeffs = P.coefficients()
Q = pol_isotop(n, k)
inds = [Pol.varname_key(str(u))[1] for u in Q.variables()]
p = 0
for mon in enumerate(P.monomials()):
m = Pol(mon[1])
p += coeffs[mon[0]]*Q.subs({x[i]:rule_sub(i, m) for i in inds})
return p
# Philip Turecek, Jun 17 2023
CROSSREFS
KEYWORD
nonn
AUTHOR
Christian G. Bower, Jan 16 2004
EXTENSIONS
a(2) corrected by Philip Turecek, Jun 13 2023
STATUS
approved