A258567 a(1) = 1; thereafter a(n) = smallest prime factor of the powerful number A001694(n).

1, 2, 2, 3, 2, 5, 3, 2, 2, 7, 2, 2, 3, 2, 2, 11, 5, 2, 2, 13, 2, 2, 2, 3, 3, 2, 2, 17, 2, 7, 19, 2, 2, 2, 3, 2, 2, 2, 23, 2, 5, 2, 3, 2, 3, 2, 2, 29, 2, 2, 31, 2, 2, 2, 2, 3, 3, 2, 2, 5, 2, 3, 11, 2, 37, 2, 2, 3, 2, 2, 41, 2, 2, 2, 43, 2, 2, 2, 3, 2, 2, 3

Offset: 1

Author

Reinhard Zumkeller, Jun 06 2015

Keywords

• nonn

Links

• Reinhard Zumkeller, Table of n, a(n) for n = 1..10000

Crossrefs

Cf. A001694, A020639, A258599, A258568, A258569, A258570, A258571.

Programs

Haskell

a258567 = a020639 . a001694

Mathematica

Table[If[Min[(f = FactorInteger[n])[[;; , 2]]] > 1 || n == 1, f[[1, 1]], Nothing], {n, 1, 3000}] (* Amiram Eldar, Jan 30 2023 *)

Python

from math import isqrt
from sympy import mobius, integer_nthroot, primefactors
def A258567(n):
    def squarefreepi(n):
        return int(sum(mobius(k)*(n//k**2) for k in range(1, isqrt(n)+1)))
    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):
        c, l = n+x, 0
        j = isqrt(x)
        while j>1:
            k2 = integer_nthroot(x//j**2, 3)[0]+1
            w = squarefreepi(k2-1)
            c -= j*(w-l)
            l, j = w, isqrt(x//k2**3)
        c -= squarefreepi(integer_nthroot(x, 3)[0])-l
        return c
    return min(primefactors(bisection(f,n,n)),default=1) # Chai Wah Wu, Sep 10 2024

Formula

a(n) = A020639(A001694(n)).

a(A258599(n)) = A000040(n) and a(m) != A000040(n) for m < A258599(n).

Extensions

Definition made more precise by N. J. A. Sloane, Apr 29 2024