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”).

A290154
Smallest number k such that exactly half the numbers in [1..k] are prime(n)-smooth.
5
6, 20, 42, 78, 118, 184, 248, 332, 428, 534, 654, 772, 906, 1052, 1208, 1388, 1562, 1754, 1958, 2164, 2396, 2638, 2896, 3144, 3424, 3682, 3986, 4304, 4622, 4976, 5286, 5652, 6002, 6374, 6748, 7148, 7532, 7934, 8356, 8786, 9224, 9684, 10158, 10618, 11114, 11604
OFFSET
1,1
COMMENTS
All terms are even numbers (because of the "exactly half the numbers in [1..k]" part of the definition).
LINKS
Michael S. Branicky, Table of n, a(n) for n = 1..10000 (terms 1..2000 from Robert Israel)
Michael S. Branicky, Python program for bfile
EXAMPLE
The 2-smooth numbers are 1, 2, 4, 8, 16, 32, ... (A000079, the powers of 2), so the numbers of 2-smooth numbers in the interval [1..k] for k = 2, 4, and 6 are 2, 3, and 3, respectively; thus, the smallest k at which the number of 2-smooth numbers in [1..k] is exactly k/2 is k=6, so a(1)=6.
The 3-smooth numbers are 1, 2, 3, 4, 6, 8, 9, 12, 16, 18, 24, 27, 32, ... (A003586), so there are more than k/2 3-smooth numbers in [1..k] for every positive k < 20, but exactly k/2 3-smooth numbers in [1..20], so a(2) = 20.
MAPLE
N:= 100:
mypi:= proc(n) option remember; global pmax; local k;
k:= procname(pmax);
while pmax < n do pmax:= nextprime(pmax); k:= k+1 od;
k
end proc:
pmax:= 2: mypi(2):= 1:
V:= Vector(N):
count:= 0:
loheap:=heap[new](`<`, 0): nlo:= 1:
hiheap:= heap[new](`>`, 1): nhi:= 1:
for k from 4 by 2 while count < N do
for v in [mypi(max(numtheory:-factorset(k-1))), mypi(max(numtheory:-factorset(k)))] do
if v <= heap[max](loheap) then heap[insert](v, loheap); nlo:= nlo+1;
elif v >= heap[max](hiheap) then heap[insert](v, hiheap); nhi:= nhi+1;
elif nlo <= nhi then heap[insert](v, loheap); nlo:= nlo+1;
else heap[insert](v, hiheap); nhi:= nhi+1;
fi;
od;
if nlo < nhi-1 then
t:= heap[extract](hiheap);
heap[insert](t, loheap);
nlo:= nlo+1; nhi:= nhi-1;
elif nhi < nlo-1 then
t:= heap[extract](loheap);
heap[insert](t, hiheap);
nhi:= nhi+1; nlo:= nlo-1;
fi;
for n from heap[max](loheap) to min(heap[max](hiheap)-1, N) do
if V[n] = 0 then count:= count+1; V[n]:= k;
fi
od;
od:
convert(V, list); # Robert Israel, Mar 28 2019
MATHEMATICA
smoothQ[k_, p_] := k <= p || Max[FactorInteger[k][[All, 1]]] <= p; a[n_] := For[p = Prime[n]; cnt = 0; k = 1, True, k++, If[smoothQ[k, p], cnt++]; If[cnt == k/2, Return[k]]]; Array[a, 46] (* Jean-François Alcover, Jul 22 2017 *)
PROG
(PARI) is(k, n) = {m=k; forprime(p=2, prime(n), while(m%p==0, m=m/p)); return(m==1); }
a(n) = {j=2; x=2; y=0; while(x!=y, j+=2; s=is(j, n)+is(j-1, n); x+=s; y+=2-s); j; } \\ Jinyuan Wang, Aug 03 2019
(Python) # see link for a faster version producing bfile
from sympy import factorint, prevprime, primerange, prod
def aupto(limit):
adict, pN = dict(), prevprime(limit+1)
pi = {p: i for i, p in enumerate(primerange(1, pN+1), start=1)}
smooth = {i: 0 for i in pi.values()}
watching = smooth[0] = 1 # 1 is prime(n) smooth for all n
for n in range(2, limit+1):
f = factorint(n, limit=pN)
nt = prod(p**f[p] for p in f if p <= pN)
if nt == n: smooth[pi[max(f)]] += 1
if 2*sum(smooth[i] for i in range(watching+1)) == n:
adict[watching] = n
watching += 1
return sorted(adict.values())
print(aupto(12000)) # Michael S. Branicky, Jun 20 2021
CROSSREFS
KEYWORD
nonn
AUTHOR
Jon E. Schoenfield, Jul 21 2017
STATUS
approved