A247246 Differences of consecutive Achilles numbers.

36, 92, 88, 104, 40, 68, 148, 27, 125, 64, 104, 4, 153, 27, 171, 29, 20, 196, 232, 144, 56, 312, 280, 108, 188, 199, 113, 67, 189, 72, 344, 16, 112, 232, 268, 63, 45, 392, 292, 32, 76, 8, 80, 587, 50, 147, 456, 184, 288, 488, 115, 772, 137, 36, 40, 212, 248

Offset: 1

Author

Eric Chen, Nov 28 2014

Keywords

• nonn
• easy

Comments

29 is the first prime in this sequence, and it equals 1352 - 1323. Clearly, if the difference is prime, then these two Achilles numbers must be relatively prime, so primes appear in this sequence rarely. However, are there infinitely many n such that a(n) is prime?

The number 1 can also appear in this sequence, because it equals 5425069448 - 5425069447 = (2^3 * 26041^2) - (7^3 * 41^2 * 97^2). Does every natural number appear in this sequence? If so, do they appear infinitely often?

Links

• Amiram Eldar, Table of n, a(n) for n = 1..10000
• Carlos Rivera, Problem 53, The Prime Puzzles and Problems Connection.
• Eric Weisstein's World of Mathematics, Achilles number.
• Wikipedia, Achilles number.

Crossrefs

Cf. A052486, A076446, A053289.

Programs

Maple

f:= proc(n) local E; E:= map(t->t[2], ifactors(n)[2]); min(E)>1 and igcd(op(E))=1 end proc:
Achilles:= select(f, [$1..10^5]):
seq(Achilles[i+1]-Achilles[i],i=1..nops(Achilles)-1); # Robert Israel, Dec 13 2014

Mathematica

achillesQ[n_] := With[{ee = FactorInteger[n][[All, 2]]}, Min[ee] > 1 && GCD @@ ee == 1];
Select[Range[10^4], achillesQ] // Differences (* Jean-François Alcover, Sep 26 2020 *)

PARI

isA052486(n) = { n>1 & vecmin(factor(n)[, 2])>1 & !ispower(n); }
lista(nn) = {v = select(n->isA052486(n), vector(nn, i, i)); vector(#v-1, n, v[n+1] - v[n]);} \\ Michel Marcus, Nov 29 2014

Python

from math import isqrt
from sympy import mobius, integer_nthroot
def A247246(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+1, 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+sum(mobius(k)*(integer_nthroot(x, k)[0]-1) for k in range(2, x.bit_length()))
        return c
    return -(a:=bisection(f,n,n))+bisection(lambda x:f(x)+1,a,a) # Chai Wah Wu, Sep 10 2024

Formula

a(n) = A052486(n+1) - A052486(n).