OFFSET
1,4
LINKS
Andrew Howroyd, Table of n, a(n) for n = 1..500
FORMULA
a(prime(n)^k) = A188445(n, k). - Andrew Howroyd, Dec 17 2018
EXAMPLE
The a(24) = 16 sets of sets with multiset union {1,1,2,3,4}:
{{1},{1,2,3,4}}
{{1,2},{1,3,4}}
{{1,3},{1,2,4}}
{{1,4},{1,2,3}}
{{1},{2},{1,3,4}}
{{1},{3},{1,2,4}}
{{1},{4},{1,2,3}}
{{1},{1,2},{3,4}}
{{1},{1,3},{2,4}}
{{1},{1,4},{2,3}}
{{2},{1,3},{1,4}}
{{3},{1,2},{1,4}}
{{4},{1,2},{1,3}}
{{1},{2},{3},{1,4}}
{{1},{2},{4},{1,3}}
{{1},{3},{4},{1,2}}
MATHEMATICA
nrmptn[n_]:=Join@@MapIndexed[Table[#2[[1]], {#1}]&, If[n==1, {}, Flatten[Cases[FactorInteger[n]//Reverse, {p_, k_}:>Table[PrimePi[p], {k}]]]]];
sqfacs[n_]:=If[n<=1, {{}}, Join@@Table[Map[Prepend[#, d]&, Select[sqfacs[n/d], Min@@#>d&]], {d, Select[Rest[Divisors[n]], SquareFreeQ]}]];
Table[Length[sqfacs[Times@@Prime/@nrmptn[n]]], {n, 90}]
PROG
(PARI)
permcount(v) = {my(m=1, s=0, k=0, t); for(i=1, #v, t=v[i]; k=if(i>1&&t==v[i-1], k+1, 1); m*=t*k; s+=t); s!/m}
sig(n)={my(f=factor(n)); concat(vector(#f~, i, vector(f[i, 2], j, primepi(f[i, 1]))))}
count(sig)={my(r=0, A=O(x*x^vecmax(sig))); for(n=1, vecsum(sig)+1, my(s=0); forpart(p=n, my(q=prod(i=1, #p, 1 + x^p[i] + A)); s+=prod(i=1, #sig, polcoef(q, sig[i]))*(-1)^#p*permcount(p)); r+=(-1)^n*s/n!); r/2}
a(n)={if(n==1, 1, my(s=sig(n)); if(#s==1, s[1]==1, count(sig(n))))} \\ Andrew Howroyd, Dec 18 2018
CROSSREFS
KEYWORD
nonn
AUTHOR
Gus Wiseman, Aug 24 2018
STATUS
approved