OFFSET
1,1
LINKS
Clark Kimberling, Table of n, a(n) for n = 1..1000
EXAMPLE
A006881 = (6, 10, 14, 15, 21, 22, 26, 33, 34, 35, 38, ... ), the increasing sequence of all products of distinct primes. The first 4 factorizations are 2*3, 2*5, 2*7, 3*5, so that (a(1), a(2), a(3), a(4)) = (2,3,4,3).
MATHEMATICA
mx = 350; t = Sort@Flatten@Table[Prime[n]*Prime[m], {n, Log[2, mx/3]}, {m, n + 1, PrimePi[mx/Prime[n]]}]; (* A006881, Robert G. Wilson v, Feb 07 2012 *)
u = Table[FactorInteger[t[[k]]][[1]], {k, 1, Length[t]}];
u1 = Table[u[[k]][[1]], {k, 1, Length[t]}] (* A096916 *)
PrimePi[u1] (* A270650 *)
v = Table[FactorInteger[t[[k]]][[2]], {k, 1, Length[t]}];
v1 = Table[v[[k]][[1]], {k, 1, Length[t]}] (* A070647 *)
PrimePi[v1] (* A270652 *)
d = v1 - u1 (* A176881 *)
Map[PrimePi[FactorInteger[#][[-1, 1]]] &, Select[Range@ 240, And[SquareFreeQ@ #, PrimeOmega@ # == 2] &]] (* Michael De Vlieger, Apr 25 2016 *)
PROG
(Python)
from math import isqrt
from sympy import primepi, primerange, primefactors
def A270652(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 int(n+x+(t:=primepi(s:=isqrt(x)))+(t*(t-1)>>1)-sum(primepi(x//k) for k in primerange(1, s+1)))
return primepi(max(primefactors(bisection(f, n, n)))) # Chai Wah Wu, Oct 23 2024
CROSSREFS
KEYWORD
nonn,easy,changed
AUTHOR
Clark Kimberling, Apr 25 2016
STATUS
approved