OFFSET
1,4
LINKS
Antti Karttunen, Table of n, a(n) for n = 1..10000
FORMULA
a(n) = A000005(n) * A067029(n)/(1+A067029(n)) = d(n) * e_n/(e_n + 1), where d(n) is the number of positive divisors of n and e_n is the exponent of the smallest prime to divide n in the prime factorization of n.
a(p) = 1 iff p is prime. - Bernard Schott, May 06 2020
a(n) = A000005(n/p) where p is the smallest prime dividing n. - David A. Corneth, May 06 2020
EXAMPLE
The divisors of 12 which are themselves divisible by 2 (the smallest prime dividing 12) are 2, 4, 6 and 12. So the 12th term is 4.
MATHEMATICA
a[1] = 0; a[n_] := DivisorSigma[0, n] * (e = FactorInteger[n][[1, 2]])/(e + 1); Array[a, 100] (* Amiram Eldar, May 06 2020 *)
PROG
(Scheme) (define (A069157 n) (let ((e_n (A067029 n))) (* (/ e_n (+ 1 e_n)) (A000005 n)))) ;; (After the formula given by the author of the sequence) - Antti Karttunen, Aug 12 2017
(Python)
from sympy import divisor_count, factorint
def a067029(n): return 0 if n==1 else next(iter(factorint(n).values()))
def a(n): return divisor_count(n)*a067029(n)//(1 + a067029(n))
print([a(n) for n in range(1, 51)]) # Indranil Ghosh, Aug 12 2017
(PARI) a(n) = if (n==1, 0, my(p=vecmin(factor(n)[, 1])); sumdiv(n, d, ((d % p) == 0))); \\ Michel Marcus, May 06 2020
CROSSREFS
KEYWORD
nonn,easy
AUTHOR
Leroy Quet, Apr 08 2002
STATUS
approved