This is a front-end for the Online Encyclopedia of Integer Sequences, made by Christian Perfect. The idea is to provide OEIS entries in non-ancient HTML, and then to think about how they're presented visually. The source code is on GitHub.
import Data.Set (Set, singleton, insert, deleteFindMin) smooth :: Set Integer -> [Integer] smooth s = x : smooth (insert (3*x) $ insert (2*x) s') where (x, s') = deleteFindMin s a003586_list = smooth (singleton 1) a003586 n = a003586_list !! (n-1) -- Reinhard Zumkeller, Dec 16 2010
[n: n in [1..4000] | PrimeDivisors(n) subset [2,3]]; // Bruno Berselli, Sep 24 2012
A003586 := proc(n) option remember; if n = 1 then 1; else for a from procname(n-1)+1 do numtheory[factorset](a) minus {2,3} ; if % = {} then return a; end if; end do: end if; end proc: # R. J. Mathar, Feb 28 2011 with(numtheory): for i from 1 to 23328 do if(i/phi(i)=3)then print(i/6) fi od; # Gary Detlefs, Jun 28 2011
a[1] = 1; j = 1; k = 1; n = 100; For[k = 2, k <= n, k++, If[2*a[k - j] < 3^j, a[k] = 2*a[k - j], {a[k] = 3^j, j++}]]; Table[a[i], {i, 1, n}] (* Hai He (hai(AT)mathteach.net) and Gilbert Traub, Dec 28 2004 *)
aa = {}; Do[If[EulerPhi[6 n] == 2 n, AppendTo[aa, n]], {n, 1, 1000}]; aa (* Artur Jasinski, Nov 05 2008 *)
fQ[n_] := Union[ MemberQ[{1, 5}, # ] & /@ Union@ Mod[ Rest@ Divisors@ n, 6]] == {False}; fQ[1] = True; Select[ Range@ 4000, fQ] (* Robert G. Wilson v, Oct 26 2010 *)
powerOfTwo = 12; Select[Nest[Union@Join[#, 2*#, 3*#] &, {1}, powerOfTwo-1], # < 2^powerOfTwo &] (* Robert G. Wilson v and T. D. Noe, Mar 03 2011 *)
fQ[n_] := n == 3 EulerPhi@ n; Select[6 Range@ 4000, fQ]/6 (* Robert G. Wilson v, Jul 08 2011 *)
mx = 4000; Sort@ Flatten@ Table[2^i*3^j, {i, 0, Log[2, mx]}, {j, 0, Log[3, mx/2^i]}] (* Robert G. Wilson v, Aug 17 2012 *)
f[n_] := Block[{p2, p3 = 3^Range[0, Floor@ Log[3, n] + 1]}, p2 = 2^Floor[Log[2, n/p3] + 1]; Min[ Select[ p2*p3, IntegerQ]]]; NestList[f, 1, 54] (* Robert G. Wilson v, Aug 22 2012 *)
Select[Range@4000, Last@Map[First, FactorInteger@#] <= 3 &] (* Vincenzo Librandi, Aug 25 2016 *)
Select[Range[4000],Max[FactorInteger[#][[All,1]]]<4&] (* Harvey P. Dale, Jan 11 2017 *)
test(n)=for(p=2,3, while(n%p==0, n/=p)); n==1; for(n=1,4000,if(test(n),print1(n",")))
list(lim)=my(v=List(),N);for(n=0,log(lim\1+.5)\log(3),N=3^n;while(N<=lim,listput(v,N);N<<=1));vecsort(Vec(v)) \\ Charles R Greathouse IV, Jun 28 2011
is_A003586(n)=n<5||vecmax(factor(n,5)[, 1])<5 \\ M. F. Hasler, Jan 16 2015
list(lim)=my(v=List(), N); for(n=0, logint(lim\=1,3), N=3^n; while(N<=lim, listput(v, N); N<<=1)); Set(v) \\ Charles R Greathouse IV, Jan 10 2018
from itertools import count, takewhile
def aupto(lim):
pows2 = list(takewhile(lambda x: xMichael S. Branicky, Jul 08 2022
from sympy import integer_log def A003586(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 f(x): return n+x-sum((x//3**i).bit_length() for i in range(integer_log(x,3)[0]+1)) return bisection(f,n,n) # Chai Wah Wu, Sep 15 2024
# faster for initial segment of sequence import heapq from itertools import islice def A003586gen(): # generator of terms v, oldv, h, psmooth_primes, = 1, 0, [1], [2, 3] 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(A003586gen(), 65))) # Michael S. Branicky, Sep 17 2024 (C++) // Returns A003586 <= threshold without approximations nor sorting #includestd::forward_list A003586(const int threshold) { std::forward_list sequence; auto start_it = sequence.before_begin(); for (int i = 1; i <= threshold; i *= 2) { for (int inc = 1; std::next(start_it) != sequence.end() && inc <= i; inc *= 3) ++start_it; auto it = start_it; for (int j = 1; i * j <= threshold; j *= 3) { sequence.emplace_after(it, i * j); for (int inc = 1; std::next(it) != sequence.end() && inc <= i; inc *= 2) ++it; } } return sequence; } // Eben Gino Lester, Apr 17 2025
def isA003586(n) : return not any(d != 2 and d != 3 for d in prime_divisors(n)) @CachedFunction def A003586(n) : if n == 1 : return 1 k = A003586(n-1) + 1 while not isA003586(k) : k += 1 return k [A003586(n) for n in (1..55)] # Peter Luschny, Jul 20 2012
import Data.Set (singleton, deleteFindMin, insert)
a033845 n = a033845_list !! (n-1)
a033845_list = f (singleton (2*3)) where
f s = m : f (insert (2*m) $ insert (3*m) s') where
(m,s') = deleteFindMin s
-- Reinhard Zumkeller, Sep 13 2011
mx = 12000; Sort@ Flatten@ Table[2^i*3^j, {i, Log[2, mx]}, {j, Log[3, mx/2^i]}] (* Robert G. Wilson v, Aug 17 2012 *)
list(lim)=my(v=List(), N); for(n=0, log(lim\2)\log(3), N=6*3^n; while(N<=lim, listput(v, N); N<<=1)); vecsort(Vec(v)) \\ Charles R Greathouse IV, Jan 02 2012
from sympy import integer_log def A033845(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 f(x): return n+x-sum((x//3**i).bit_length() for i in range(integer_log(x,3)[0]+1)) return 6*bisection(f,n,n) # Chai Wah Wu, Sep 15 2024
# faster for initial segment of sequence import heapq from itertools import islice def A033845gen(): # generator of terms v, oldv, h, psmooth_primes, = 1, 0, [1], [2, 3] while True: v = heapq.heappop(h) if v != oldv: yield 6*v oldv = v for p in psmooth_primes: heapq.heappush(h, v*p) print(list(islice(A033845gen(), 50))) # Michael S. Branicky, Sep 18 2024
a022329 n = a022329_list !! (n-1) -- Where a022329_list is defined in A022328. -- Reinhard Zumkeller, Nov 19 2015, May 16 2015
IntegerExponent[Select[Range[10^5], # == 2^IntegerExponent[#, 2] * 3^IntegerExponent[#, 3] &], 3] (* Amiram Eldar, Apr 15 2024 *)
from sympy import integer_log def A022329(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 f(x): return n+x-sum((x//3**i).bit_length() for i in range(integer_log(x,3)[0]+1)) return integer_log((m:=bisection(f,n,n))>>(~m&m-1).bit_length(),3)[0] # Chai Wah Wu, Sep 15 2024
a(10) = 2 because A033845(10) = 108 = 2^2*3^3. a(100) = 2 because A033845(100) = 59872 = 2^8*3^7. a(1000) = 56 because A033845(1000) = 216172782113783808 = 2^56*3^1.
mx = 1000000; t = Select[Sort[Flatten[Table[2^i 3^j, {i, Log[2, mx]}, {j, Log[3, mx]}]]], # <= mx &]; Table[FactorInteger[i][[1, 2]], {i, t}] (* T. D. Noe, Aug 31 2012 *)
from sympy import integer_log def A191475(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 f(x): return n+x-sum((x//3**i).bit_length() for i in range(integer_log(x,3)[0]+1)) return 1+(~(m:=bisection(f,n,n))&m-1).bit_length() # Chai Wah Wu, Sep 15 2024
a(10) = 5 because A033845(10) = 108 = 2^2 + 3^3 => A191475(10) = 2, A191476(10) = 3.
Comments