A097621 In canonical prime factorization of n replace p^e with its index in A000961.

1, 2, 3, 4, 5, 6, 6, 7, 8, 10, 9, 12, 10, 12, 15, 11, 12, 16, 13, 20, 18, 18, 14, 21, 15, 20, 16, 24, 17, 30, 18, 19, 27, 24, 30, 32, 20, 26, 30, 35, 21, 36, 22, 36, 40, 28, 23, 33, 24, 30, 36, 40, 25, 32, 45, 42, 39, 34, 26, 60, 27, 36, 48, 28, 50, 54, 29, 48, 42, 60, 30, 56, 31

Offset: 1

Author

Reinhard Zumkeller, Aug 17 2004

Keywords

• nonn
• mult

Comments

The definition of the sequence has been corrected, given that it uses A095874, the indices in the list A000961 of "powers of primes" starting with A000961(1) = 1, rather than A322981, index of p^e in the list of prime powers A246655, as written in the original definition. See A333235 for the variant of this sequence which uses A322981 and A246655 instead, maybe the more natural choice given that the factorization of integers consists of prime powers > 1. - M. F. Hasler, Jun 15 2021

Examples

			n=600 = 2^3 * 3 * 5^2 -> A095874(8)*A095874(3)*A095874(25) = 7 * 3 * 15 = 315.
		

Links

• Ivan Neretin, Table of n, a(n) for n = 1..10000

Crossrefs

Cf. A000961 (powers of primes), A246655 (prime powers), A003963, A018266, A095874 (index of n = p^e in A000961).

Cf. A097622, A097623, A097624.

Cf. A322981 (index of n = p^e in A246655), A333235 (variant of this sequence).

Programs

Maple

N:= 1000: # to get a(1) to a(N)
Primes:= select(isprime,[2, seq(2*i+1,i=1..(N-1)/2)]):
PP:= sort([1,seq(seq(p^k, k=1..floor(log[p](N))),p=Primes)]):
for n from 1 to nops(PP) do B[PP[n]]:= n od:
seq(mul(B[f[1]^f[2]],f=ifactors(n)[2]),n=1..N); # Robert Israel, Sep 02 2015

Mathematica

pp = Select[Range@100, Length[FactorInteger[#]] == 1 &]; a = Table[Times @@ (Position[pp, #][[1, 1]] & /@ (#[[1]]^#[[2]] & /@ FactorInteger[n])), {n, 73}] (* Ivan Neretin, Sep 02 2015 *)

PARI

f(n) = if(isprimepower(n), sum(i=1, logint(n, 2), primepi(sqrtnint(n, i)))+1, n==1); \\ A095874
a(n) = my(fr=factor(n)); for (k=1, #fr~, fr[k,1] = f(fr[k,1]^fr[k,2]); fr[k,2] = 1); factorback(fr); \\ Michel Marcus, May 29 2021
A097621(n)=vecprod([A095874(f[1]^f[2])|f<-factor(n)~]) \\ M. F. Hasler, Jun 15 2021

Python

from math import prod
from sympy import primepi, integer_nthroot, factorint
def A097621(n): return prod(1+int(primepi(m:=p**e)+sum(primepi(integer_nthroot(m,k)[0]) for k in range(2,m.bit_length()))) for p,e in factorint(n).items()) # Chai Wah Wu, Jan 19 2025

Formula

Multiplicative with: p^e -> A095874(p^e), p prime.

a(A000961(n)) = n; a(a(n)) = A097622(n); a(a(a(n))) = A097623(n);

a(n) <= n; a(n) = n iff 60 mod n = 0: a(A018266(n)) = A018266(n);

a(A097624(n)) = n and a(m) <> n for n < A097624(n).

Extensions

Definition corrected by M. F. Hasler, Jun 16 2021

Example corrected by Ray Chandler, Jun 30 2021