A137484 - OEIS

%I #24 Feb 22 2025 01:22:57

%S 576,1600,2916,3136,7744,10816,18225,18496,23104,33856,35721,53824,

%T 61504,62500,87616,88209,107584,118336,123201,140625,141376,179776,

%U 210681,222784,238144,263169,287296,322624,341056,385641,399424,440896

%N Numbers with 21 divisors.

%C Maple implementation: see A030513.

%C Numbers of the form p^20 or p^2*q^6 (A189990) where p and q are distinct primes. - _R. J. Mathar_, Mar 01 2010

%H T. D. Noe, <a href="/A137484/b137484.txt">Table of n, a(n) for n = 1..1000</a>

%F A000005(a(n)) = 21.

%F Sum_{n>=1} 1/a(n) = P(2)*P(6) - P(8) + P(20) = 0.00365945..., where P is the prime zeta function. - _Amiram Eldar_, Jul 03 2022

%t Select[Range[450000],DivisorSigma[0,#]==21&] (* _Vladimir Joseph Stephan Orlovsky_, May 03 2011 *)

%o (PARI) is(n)=numdiv(n)==21 \\ _Charles R Greathouse IV_, Jun 19 2016

%o (Python)

%o from math import isqrt

%o from sympy import primepi, primerange, integer_nthroot

%o def A137484(n):

%o     def bisection(f,kmin=0,kmax=1):

%o         while f(kmax) > kmax: kmax <<= 1

%o         kmin = kmax >> 1

%o         while kmax-kmin > 1:

%o             kmid = kmax+kmin>>1

%o             if f(kmid) <= kmid:

%o                 kmax = kmid

%o             else:

%o                 kmin = kmid

%o         return kmax

%o     def f(x): return n+x-sum(primepi(isqrt(x//p**6)) for p in primerange(integer_nthroot(x,6)[0]+1))+primepi(integer_nthroot(x,8)[0])-primepi(integer_nthroot(x,20)[0])

%o     return bisection(f,n,n) # _Chai Wah Wu_, Feb 21 2025

%Y Cf. A000005, A030513, A030638 (20 divisors), A137485 (22 divisors), A189990.

%K nonn,changed

%O 1,1

%A _R. J. Mathar_, Apr 22 2008