A188585 Moebius inversion of sequence A000688, the number of factorizations of n into prime powers greater than 1.

1, 0, 0, 1, 0, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 0, 0, 2, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 4, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0

Offset: 1

Author

Marc Bogaerts, Apr 04 2011

Keywords

• nonn
• mult
• easy

Comments

Dirichlet convolution product of A000688 with the Moebius function.

It appears that a(n) is nonzero for n in A001694, the powerful numbers. - T. D. Noe, Apr 06 2011 [This is correct: a(n) > 0 if and only if n is in A001694. - Amiram Eldar, Jun 10 2025]

There is a similar sequence defined by b(n) = Product_{i} floor(e(i)/2) where n = Product_{p} p(i)^e(i) is the usual prime factorization, which differs from a(n) at n = 64, 128, 256, 512, 576, 729,.... - R. J. Mathar, Sep 18 2012 [This sequence is A365550. - Amiram Eldar, Jun 10 2025]

The number of unordered factorizations of n into 1 and prime powers p^e where p is prime and e >= 2 (A025475). - Amiram Eldar, Jun 10 2025

Links

• Charles R Greathouse IV, Table of n, a(n) for n = 1..10000

Crossrefs

Cf. A000688, A001694, A002865, A008683, A025475, A365550.

Programs

GAP

mtrf:=function ( f, x )     # the Moebius inversion formula
    local  d;
    d := DivisorsInt( x );
    return Sum( d, function ( i )
            return f( i ) * MoebiusMu( (x / i) );
        end );
end;
nra:=function ( x )         # the number of Abelian groups of order x
    local  pp, ll;
    pp := PrimePowersInt( x );
    ll := [ 1 .. Size( pp ) / 2 ];
    return Product( List( 2 * ll, function ( i )
              return NrPartitions( pp[i] );
          end ) );
end;
a:=function ( n )
    return mtrf( nra, n );
end;

Maple

with(numtheory): with(combinat):
a:= n-> add(mobius(n/d) *mul(numbpart(i[2]),
        i=ifactors(d)[2]), d=divisors(n)):
seq(a(n), n=1..110);  # Alois P. Heinz, Apr 07 2011

Mathematica

MobiusTransform[a_List] := Module[{n = Length[a], b}, b = Table[0, {i, n}];Do[b[[i]] = Plus @@ (MoebiusMu[i/Divisors[i]] a[[Divisors[i]]]), {i, n}]; b]; A688[n_] := Times @@ PartitionsP /@ Last /@ FactorInteger@n; MobiusTransform[Array[A688, 100]] (* T. D. Noe, Apr 06 2011 *)
f[p_, e_] := PartitionsP[e] - PartitionsP[e-1]; a[1] = 1; a[n_] := Times @@ f @@@ FactorInteger[n]; Array[a, 100] (* Amiram Eldar, Jun 10 2025 *)

PARI

a(n) = vecprod(apply(x -> numbpart(x)-numbpart(x-1), factor(n)[, 2])); \\ Amiram Eldar, Jun 10 2025

Python

from math import prod
from sympy import partition, factorint
def A188585(n): return prod(partition(e)-partition(e-1) for e in factorint(n).values()) # Chai Wah Wu, Jun 10 2025

Formula

a(n) = Sum_{d|n} A008683(n/d) * A000688(d).

Dirichlet g.f.: Product_{k>=2} zeta(k*s). - Ilya Gutkovskiy, Nov 03 2020

Sum_{k=1..n} a(k) ~ c * sqrt(n), where c = Product_{k>=3} zeta(k/2) = 10.0301441966843566206076085895839492473559217336... - Vaclav Kotesovec, Apr 22 2025

Multiplicative with a(p^e) = A002865(e). - Amiram Eldar, Jun 10 2025