%I #16 May 19 2019 23:36:20
%S 8,9,18,24,25,32,36,40,45,49,50,56,63,64,75,81,88,90,96,98,99,100,104,
%T 117,120,121,125,126,128,136,144,147,150,152,153,160,162,168,169,171,
%U 175,180,184,192,196,198,200,207,216,224,225,232,234,242,243,245,248
%N Integers n such that the orbit n, f(n), f(f(n)), ... is eventually periodic with period 2, where f(n) = product(a(k)^p(k)) when n has the prime factorization n = product(p(k)^a(k)).
%C It is proved in the reference that for every positive integer n the orbit n, f(n), f(f(n)), ... is eventually periodic with period 1 or 2.
%C Includes all numbers whose prime exponents are distinct primes. If n is in this sequence and k is a squarefree number such that (k,n) = 1, then k*n is in this sequence. - _Charlie Neder_, May 16 2019
%H Paolo P. Lava, <a href="/A159836/b159836.txt">Table of n, a(n) for n = 1..5000</a>
%H M. Farrokhi, <a href="https://www.jstor.org/stable/40391138">The Prime Exponentiation of an Integer: Problem 11315</a>, Amer. Math. Monthly, 116 (2009), 470.
%p with(numtheory); P:=proc(q) local a0f,a1,a1f,a2,a2f,a3,a3f,a4,a4f,k,n;
%p for n from 1 to q do a0:=1;a1:=1;a2:=2;a3:=3;a4:=n;
%p while not (a1=a3 and a2=a4) do a0f:=ifactors(a4)[2];
%p a1:=mul(a0f[k][2]^a0f[k][1],k=1..nops(a0f)); a1f:=ifactors(a1)[2];
%p a2:=mul(a1f[k][2]^a1f[k][1],k=1..nops(a1f)); a2f:=ifactors(a2)[2];
%p a3:=mul(a2f[k][2]^a2f[k][1],k=1..nops(a2f)); a3f:=ifactors(a3)[2];
%p a4:=mul(a3f[k][2]^a3f[k][1],k=1..nops(a3f)); od;
%p if a1<>a2 then print(n); fi; od; end: P(10^6); # _Paolo P. Lava_, Oct 24 2013
%t f[n_] := Module[{f = Transpose[FactorInteger[n]]}, Times @@ (f[[2]]^f[[1]])]; Select[Range[300], (x = NestWhileList[f, #, UnsameQ, All]; x[[-2]] != x[[-1]]) &] (* _T. D. Noe_, Oct 24 2013 *)
%Y Cf. A008477, A008478.
%K nonn
%O 1,1
%A _John W. Layman_, Apr 23 2009