A376120 Refactorable numbers that are perfect powers.

1, 8, 9, 36, 128, 225, 441, 625, 1089, 1521, 2025, 2601, 3249, 3600, 4761, 5625, 6561, 7569, 7776, 8100, 8649, 10000, 12321, 15129, 16641, 19881, 21952, 22500, 25281, 26244, 28224, 31329, 32400, 32768, 33489, 35721, 40401, 45369, 47961, 50625, 56169, 62001, 64000, 71289, 84681, 90000

Offset: 1

Author

Waldemar Puszkarz, Sep 11 2024

Keywords

• nonn

Comments

Intersection of A001597 and A033950.

Examples

			8 is a perfect power, as 8=2^3, and it is also a refactorable numbers, being divisible by its number of divisors (4).
		

Crossrefs

Cf. A001597 (perfect powers), A033950 (refactorable numbers).

Programs

Mathematica

Join[{1}, Select[Range[10^5], Divisible[#, DivisorSigma[0,#]]&&GCD@@FactorInteger[#][[All, 2]]>1&]]

PARI

ok(n) = n==1 || (n%numdiv(n)==0&&ispower(n))

Python

from itertools import count, islice
from math import prod
from sympy import mobius, integer_nthroot, factorint
def A376120_gen(): # generator of terms
    def bisection(f,kmin=0,kmax=1):
        while f(kmax) > kmax: kmax <<= 1
        while kmax-kmin > 1:
            kmid = kmax+kmin>>1
            if f(kmid) <= kmid:
                kmax = kmid
            else:
                kmin = kmid
        return kmax
    def f(x): return int(x-1+sum(mobius(k)*(integer_nthroot(x,k)[0]-1) for k in range(2,x.bit_length())))
    m = 1
    for n in count(1):
        m = bisection(lambda x:f(x)+n,m,m)
        if not m%prod(e+1 for e in factorint(m).values()): yield m
A376120_list = list(islice(A376120_gen(),40)) # Chai Wah Wu, Oct 04 2024