with(combinat); with(numtheory); pet_cycleind_symm := proc(n) option remember; if n=0 then return 1; fi; expand(1/n*add(a[l]*pet_cycleind_symm(n-l), l=1..n)); end; TI_mu := proc(pv) local p, n, beta, ind, res; if pv = 1 then return 2 fi; p := op(2, ifactors(pv)); n := nops(p); res := 0; for ind from 2^n to 2^(n+1)-1 do beta := convert(ind, base, 2); res := res + (-1)^add(beta[q], q=1..n) *2^mul(p[q][1]^(p[q][2]-beta[q]), q=1..n); od; res/pv; end; TI_rcyc := r -> mul(a[p]^TI_mu(p), p in divisors(r)); CART_prod := proc(t1, t2) local v1, v2, l1, l2, res; res := 1; for v1 in indets(t1) do l1 := op(1, v1); for v2 in indets(t2) do l2 := op(1, v2); res := res * a[lcm(l1, l2)]^ (gcd(l1, l2) *degree(t1, v1) *degree(t2, v2)); od; od; res; end; CART_pow := proc(t, q) local res, p; res := t; for p to q-1 do res := CART_prod(res, t); od; res; end; pet_cycleind_hypergraph := proc(n) option remember; local term, v, contr, edgidx, edgterm; if n = 0 then return a[1] fi; if n = 1 then return a[1]^2 fi; edgidx := 0; for term in pet_cycleind_symm(n) do edgterm := 1; for v in indets(term) do contr := CART_pow(TI_rcyc(op(1, v)), degree(term, v)); if type(edgterm, `integer`) then edgterm := contr; else edgterm := CART_prod(edgterm, contr); fi; od; edgidx := edgidx + lcoeff(term) * edgterm; od; edgidx; end;