%I #27 Aug 14 2024 10:04:15
%S 4,9,16,25,32,36,49,81,100,121,128,144,169,196,225,243,256,289,324,
%T 361,400,441,484,529,576,625,676,784,841,900,961,1024,1089,1156,1225,
%U 1296,1369,1444,1521,1600,1681,1764,1849,1936,2025,2048,2116,2187,2209,2304,2401,2500
%N Noncube perfect powers.
%C This was the original definition of A239870. However, the true definition of that sequence seems to be slightly different.
%H Hugo Pfoertner, <a href="/A340585/b340585.txt">Table of n, a(n) for n = 1..10000</a>
%F Sum_{n>=1} 1/a(n) = 1 - zeta(3) + Sum_{k>=2} mu(k)*(1-zeta(k)) = 1 - A002117 + A072102 = 0.6724074652... - _Amiram Eldar_, Jan 12 2021
%p filter:= proc(n) local g;
%p g:= igcd(op(ifactors(n)[2][..,2]));
%p g > 1 and (g mod 3 <> 0)
%p end proc:
%p select(filter, [$1..10000]); # _Robert Israel_, Jan 12 2021
%t Select[Range[2, 2500], (g = GCD @@ FactorInteger[#][[;; , 2]]) > 1 && !Divisible[g, 3] &] (* _Amiram Eldar_, Jan 12 2021 *)
%o (PARI) for(n=2,2500,if( ispower(n) % 3, print1(n,", ")))
%o (Python)
%o from math import isqrt
%o from sympy import mobius, integer_nthroot
%o def A340585(n):
%o def f(x): return int(n+x-isqrt(x)+sum(mobius(k)*(integer_nthroot(x,k)[0]-1) for k in range(5,x.bit_length())))
%o kmin, kmax = 1,2
%o while f(kmax) >= kmax:
%o kmax <<= 1
%o while True:
%o kmid = kmax+kmin>>1
%o if f(kmid) < kmid:
%o kmax = kmid
%o else:
%o kmin = kmid
%o if kmax-kmin <= 1:
%o break
%o return kmax # _Chai Wah Wu_, Aug 14 2024
%Y Cf. A000578, A001597, A002117, A007412, A072102, A076467, A097054, A239728.
%K nonn
%O 1,1
%A _Hugo Pfoertner_, Jan 12 2021