OFFSET
1,1
COMMENTS
Terms of A216427 that are not 5th powers of squarefree numbers (A113850) and not 10th powers of primes (A030629). - Amiram Eldar, Feb 07 2023
LINKS
Charles R Greathouse IV, Table of n, a(n) for n = 1..10000
FORMULA
Sum_{n>=1} 1/a(n) = 2 + ((zeta(2)-1)*(zeta(3)-1)-1)/zeta(6) - zeta(5)/zeta(10) - P(6) - P(10) = 0.09117811499514578262..., where P(s) is the prime zeta function. - Amiram Eldar, Feb 07 2023
MATHEMATICA
With[{upto=4000}, Select[Union[Flatten[{#[[1]]^2 #[[2]]^3, #[[2]]^2 #[[1]]^3}& /@ Subsets[Range[2, Surd[upto, 2]], {2}]]], #<=upto&]](* Harvey P. Dale, Jan 04 2014 *)
pMx = 25; mx = 2^3 pMx^2; t = Flatten[Table[If[a != b, a^2 b^3, 0], {a, 2, mx^(1/2)}, {b, 2, mx^(1/3)}]]; Union[Select[t, 0 < # <= mx &]] (* T. D. Noe, Jan 02 2014 *)
PROG
(PARI) list(lim)=my(v=List()); for(b=2, sqrtnint(lim\4, 3), for(a=2, sqrtint(lim\b^3), if(a!=b, listput(v, a^2*b^3)))); Set(v) \\ Charles R Greathouse IV, Jan 02 2014
(Python)
from math import isqrt
from sympy import integer_nthroot, mobius, primepi
def A216426(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):
j, b, a, d = isqrt(x), integer_nthroot(x, 6)[0], integer_nthroot(x, 5)[0], integer_nthroot(x, 10)[0]
l, c = 0, n+x-2+primepi(b)+sum(mobius(k)*(j//k**3) for k in range(d+1, b+1))+primepi(d)+sum(mobius(k)*(a//k**2+j//k**3) for k in range(1, d+1))
while j>1:
k2 = integer_nthroot(x//j**2, 3)[0]+1
w = sum(mobius(k)*((k2-1)//k**2) for k in range(1, isqrt(k2-1)+1))
c -= j*(w-l)
l, j = w, isqrt(x//k2**3)
return c+l
return bisection(f, n, n) # Chai Wah Wu, Sep 13 2024
CROSSREFS
KEYWORD
nonn
AUTHOR
V. Raman, Sep 07 2012
EXTENSIONS
Name corrected by Charles R Greathouse IV, Jan 02 2014
STATUS
approved