cp's OEIS Frontend

A133581 (k^2)-th k-smooth number for k = prime(n).

%I A133581 #14 Feb 16 2025 08:33:06
%S A133581 8,16,54,112,396,512,1008,1155,1794,3312,3520,5488,6776,7020,8405,
%T A133581 11180,14384,14720,18241,20339,20709,24769,27094,31648,38994,41890,
%U A133581 42336,45318,45825,48852,66234,69874,76857,77441,91719,92323,100215,108376,112896,121539
%N A133581 (k^2)-th k-smooth number for k = prime(n).
%C A133581 An integer is k-smooth if it has no prime factors > k.
%H A133581 Eric Weisstein's World of Mathematics, <a href="https://mathworld.wolfram.com/SmoothNumber.html">Smooth Number</a>
%F A133581 a(n) = A001248(n)-th integer which has no prime factors > A000040(n).
%e A133581 a(1) = 8 = A000079(4).
%e A133581 a(2) = 16 = A003586(9).
%e A133581 a(3) = 54 = A051037(25).
%o A133581 (Python)
%o A133581 from sympy import integer_log, prime, prevprime
%o A133581 def A133581(n):
%o A133581     if n==1: return 8
%o A133581     def bisection(f,kmin=0,kmax=1):
%o A133581         while f(kmax) > kmax: kmax <<= 1
%o A133581         while kmax-kmin > 1:
%o A133581             kmid = kmax+kmin>>1
%o A133581             if f(kmid) <= kmid:
%o A133581                 kmax = kmid
%o A133581             else:
%o A133581                 kmin = kmid
%o A133581         return kmax
%o A133581     def g(x,m): return sum((x//3**i).bit_length() for i in range(integer_log(x,3)[0]+1)) if m==3 else sum(g(x//(m**i),prevprime(m))for i in range(integer_log(x,m)[0]+1))
%o A133581     k = prime(n)
%o A133581     def f(x): return k**2+x-g(x,k)
%o A133581     return bisection(f,k**2,k**2) # _Chai Wah Wu_, Sep 17 2024
%Y A133581 Cf. A000040, A000079, A001248, A003586, A051037, A002473, A051038.
%K A133581 nonn,less
%O A133581 1,1
%A A133581 _Jonathan Vos Post_, Dec 26 2007
%E A133581 Corrected and extended by _D. S. McNeil_, Dec 08 2010
%E A133581 a(33)-a(40) from _Chai Wah Wu_, Sep 17 2024