login

Year-end appeal: Please make a donation to the OEIS Foundation to support ongoing development and maintenance of the OEIS. We are now in our 61st year, we have over 378,000 sequences, and we’ve reached 11,000 citations (which often say “discovered thanks to the OEIS”).

A080197
13-smooth numbers: numbers whose prime divisors are all <= 13.
25
1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 18, 20, 21, 22, 24, 25, 26, 27, 28, 30, 32, 33, 35, 36, 39, 40, 42, 44, 45, 48, 49, 50, 52, 54, 55, 56, 60, 63, 64, 65, 66, 70, 72, 75, 77, 78, 80, 81, 84, 88, 90, 91, 96, 98, 99, 100, 104, 105, 108, 110, 112, 117, 120
OFFSET
1,2
COMMENTS
Numbers of the form 2^r*3^s*5^t*7^u*11^v*13^w with r, s, t, u, v, w >= 0.
LINKS
FORMULA
Sum_{n>=1} 1/a(n) = Product_{primes p <= 13} p/(p-1) = (2*3*5*7*11*13)/(1*2*4*6*10*12) = 1001/192. - Amiram Eldar, Sep 22 2020
EXAMPLE
33 = 3*11 and 39 = 3*13 are terms but 34 = 2*17 is not.
MATHEMATICA
mx = 120; Sort@ Flatten@ Table[ 2^i*3^j*5^k*7^l*11^m*13^n, {i, 0, Log[2, mx]}, {j, 0, Log[3, mx/2^i]}, {k, 0, Log[5, mx/(2^i*3^j)]}, {l, 0, Log[7, mx/(2^i*3^j*5^k)]}, {m, 0, Log[11, mx/(2^i*3^j*5^k*7^l)]}, {n, 0, Log[13, mx/(2^i*3^j*5^k*7^l*11^m)]}] (* Robert G. Wilson v, Aug 17 2012 *)
PROG
(PARI) test(n)=m=n; forprime(p=2, 13, while(m%p==0, m=m/p)); return(m==1)
for(n=1, 200, if(test(n), print1(n", ")))
(PARI) is_A080197(n, p=13)=n<=p||vecmax(factor(n, p+1)[, 1])<=p \\ M. F. Hasler, Jan 16 2015
(PARI) list(lim, p=13)=if(p==2, return(powers(2, logint(lim\1, 2)))); my(v=[], q=precprime(p-1), t=1); for(e=0, logint(lim\=1, p), v=concat(v, list(lim\t, q)*t); t*=p); Set(v) \\ Charles R Greathouse IV, Apr 16 2020
(Magma) [n: n in [1..150] | PrimeDivisors(n) subset PrimesUpTo(13)]; // Bruno Berselli, Sep 24 2012
(Python)
import heapq
from itertools import islice
from sympy import primerange
def agen(p=13): # generate all p-smooth terms
v, oldv, h, psmooth_primes, = 1, 0, [1], list(primerange(1, p+1))
while True:
v = heapq.heappop(h)
if v != oldv:
yield v
oldv = v
for p in psmooth_primes:
heapq.heappush(h, v*p)
print(list(islice(agen(), 69))) # Michael S. Branicky, Nov 20 2022
(Python)
from sympy import integer_log, prevprime
def A080197(n):
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 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))
def f(x): return n+x-g(x, 13)
return bisection(f, n, n) # Chai Wah Wu, Sep 16 2024
CROSSREFS
Cf. A000079, A080196. For p-smooth numbers with other values of p, see A003586, A051037, A002473, A051038, A080681, A080682, A080683.
Sequence in context: A102800 A072676 A235986 * A341749 A115847 A204315
KEYWORD
easy,nonn
AUTHOR
Klaus Brockhaus, Feb 10 2003
STATUS
approved