A385133 The number of integers k from 1 to n such that gcd(n, k) has an odd maximum exponent in its prime factorization.

0, 1, 1, 1, 1, 4, 1, 3, 2, 6, 1, 5, 1, 8, 7, 5, 1, 10, 1, 7, 9, 12, 1, 13, 4, 14, 7, 9, 1, 22, 1, 11, 13, 18, 11, 12, 1, 20, 15, 19, 1, 30, 1, 13, 16, 24, 1, 23, 6, 28, 19, 15, 1, 32, 15, 25, 21, 30, 1, 29, 1, 32, 20, 21, 17, 46, 1, 19, 25, 46, 1, 33, 1, 38, 32

Offset: 1

Author

Amiram Eldar, Jun 24 2025

Keywords

• nonn
• easy

Links

• Amiram Eldar, Table of n, a(n) for n = 1..10000

Crossrefs

Cf. A001221, A051903, A368714, A384655, A385132.

Programs

Mathematica

e[n_] := If[n == 1, 0, Max[FactorInteger[n][[;; , 2]]]]; a[n_] := Sum[Boole[OddQ[e[GCD[n, k]]]], {k, 1, n}]; Array[a, 100]
(* second program: *)
a[n_] := Module[{f = FactorInteger[n], p, e, emax, kmax}, p = f[[;;, 1]]; e = f[[;;, 2]]; emax = Max[e]; kmax = emax + Mod[emax, 2]; Sum[(-1)^k * Product[p[[i]]^e[[i]] - If[e[[i]] < k, 0, p[[i]]^(e[[i]]-k)], {i, 1, Length[p]}], {k, 1, kmax}]]; a[1] = 0; Array[a, 100]

PARI

q(n) = if(n == 1, 0, vecmax(factor(n)[,2]) % 2);
a(n) = sum(k = 1, n, q(gcd(n, k)));

PARI

a(n) = if(n == 1, 0, my(f = factor(n), p = f[,1], e = f[,2], emax = vecmax(e), kmax = emax + emax % 2); sum(k = 1, kmax, (-1)^k * prod(i = 1, #p, p[i]^e[i] - if(e[i] < k, 0, p[i]^(e[i]-k)))));

Formula

a(n) = Sum_{k=1..n} (A051903(gcd(n, k)) mod 2).

a(n) = n - A385132(n).

a(n) = Sum_{k=1..kmax(n)} (-1)^k * Product_{i=1..r} f(p_i^e_i, k), for n >= 2; if n = Product_{i=1..r} p_i^e_i, r = omega(n) = A001221(n), then emax(n) = max(e_i) = A051903(n), kmax(n) = emax(n)+1 if emax(n) is odd, and emax(n) otherwise, f(p^e, k) = p^e - p^(e-k) if e >= k, and f(p^e, k) = p^e if e < k.

Sum_{k=1..n} a(k) ~ c * n^2 / 2, where c = Sum_{k>=1} (-1)^(k+1)*(1-1/zeta(2*k)) = 0.32979448728905469699... .