OFFSET
1,1
COMMENTS
Number of terms < 10^n: 0, 1, 13, 60, 252, 916, 3158, 10553, 34561, 111891, 359340, 1148195, 3656246, 11616582, 36851965, ..., A118896(n) - A070428(n). - Robert G. Wilson v, Aug 11 2014
a(n) = (s(n))^2 * f(n), s(n) > 1, f(n) > 1, where s(n) is not a power of f(n), and f(n) is squarefree and gcd(s(n), f(n)) = f(n). - Daniel Forgues, Aug 11 2015
LINKS
Amiram Eldar, Table of n, a(n) for n = 1..10000 (terms 1..1000 from T. D. Noe)
Robert Israel, log-log plot of a(n).
OEIS Wiki, Achilles numbers.
Project Euler, Problem 302: Strong Achilles Numbers.
Eric Weisstein's World of Mathematics, Achilles Number.
Chai Wah Wu, Algorithms for Complementary Sequences, Integers (2025) Vol. 25, Art. No. A95. See p. 24.
FORMULA
a(n) = O(n^2). - Daniel Forgues, Aug 11 2015
a(n) = O(n^2 / log log n). - Daniel Forgues, Aug 12 2015
Sum_{n>=1} 1/a(n) = zeta(2)*zeta(3)/zeta(6) - Sum_{k>=2} mu(k)*(1-zeta(k)) - 1 = A082695 - A072102 - 1 = 0.06913206841581433836... - Amiram Eldar, Oct 14 2020
EXAMPLE
a(3)=200 because 200=2^3*5^2, both 3 and 2 are greater than 1, and the highest common factor of 3 and 2 is 1.
Factorizations of a(1) to a(20):
72 = 2^3 3^2, 108 = 2^2 3^3, 200 = 2^3 5^2, 288 = 2^5 3^2,
392 = 2^3 7^2, 432 = 2^4 3^3, 500 = 2^2 5^3, 648 = 2^3 3^4,
675 = 3^3 5^2, 800 = 2^5 5^2, 864 = 2^5 3^3, 968 = 2^3 11^2,
972 = 2^2 3^5, 1125 = 3^2 5^3, 1152 = 2^7 3^2, 1323 = 3^3 7^2,
1352 = 2^3 13^2, 1372 = 2^2 7^3, 1568 = 2^5 7^2, 1800 = 2^3 3^2 5^2.
Examples for a(n) = (s(n))^2 * f(n): (see above comment)
s(n) = 6, 6, 10, 12, 14, 12, 10, 18, 15, 20, 12, 22, 18, 15, 24, 21,
f(n) = 2, 3, 2, 2, 2, 3, 5, 2, 3, 2, 6, 2, 3, 5, 2, 3,
MAPLE
filter:= proc(n) local E; E:= map(t->t[2], ifactors(n)[2]); min(E)>1 and igcd(op(E))=1 end proc:
select(filter, [$1..10000]); # Robert Israel, Aug 11 2014
MATHEMATICA
achillesQ[n_] := Block[{ls = Last /@ FactorInteger@n}, Min@ ls > 1 == GCD @@ ls]; Select[ Range@ 5500, achillesQ@# &] (* Robert G. Wilson v, Jun 10 2010 *)
PROG
(PARI) is(n)=my(f=factor(n)[, 2]); n>9 && vecmin(f)>1 && gcd(f)==1 \\ Charles R Greathouse IV, Sep 18 2015, replacing code by M. F. Hasler, Sep 23 2010
(Python)
from math import gcd
from itertools import count, islice
from sympy import factorint
def A052486_gen(startvalue=1): # generator of terms >= startvalue
return (n for n in count(max(startvalue, 1)) if (lambda x: all(e > 1 for e in x) and gcd(*x) == 1)(factorint(n).values()))
(Python)
from math import isqrt
from sympy import mobius, integer_nthroot
def A052486(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 bisection(f, n, n) # Chai Wah Wu, Sep 10 2024
(SageMath)
def isA052486(n: int) -> bool:
factors = factor(n)
exps = [e for (_, e) in factors]
if any(e < 2 for e in exps): return False
return reduce(gcd, exps) == 1
print([n for n in range(2, 6100) if isA052486(n)]) # Peter Luschny, Dec 02 2025
CROSSREFS
KEYWORD
nonn
AUTHOR
Henry Bottomley, Mar 16 2000
EXTENSIONS
Example edited by Mac Coombe (mac.coombe(AT)gmail.com), Sep 18 2010
Name edited by M. F. Hasler, Jul 17 2019
STATUS
approved