OFFSET
1,2
COMMENTS
a(m) mod prime(n) > 0 for m < A258603(n); a(A258600(n)) = A030516(n) = prime(n)^6. - Reinhard Zumkeller, Jun 06 2015
LINKS
Alois P. Heinz, Table of n, a(n) for n = 1..10000
FORMULA
Sum_{n>=1} 1/a(n) = Product_{p prime} (1 + 1/(p^5*(p-1))) = 1.0334657852594050612296726462481884631303137561267151463866539131591664... - Amiram Eldar, Jul 09 2020
EXAMPLE
2^7*3^6 = 93312 is a member (although not of A076470).
MATHEMATICA
Join[{1}, Select[Range[900000], Min[FactorInteger[#][[All, 2]]]>5&]] (* Harvey P. Dale, Mar 03 2018 *)
PROG
(PARI) for(n=1, 560000, if(n*sumdiv(n, d, isprime(d)/d^6)==floor(n*sumdiv(n, d, isprime(d)/d^6)), print1(n, ", ")))
(Haskell)
import Data.Set (singleton, deleteFindMin, fromList, union)
a069493 n = a069493_list !! (n-1)
a069493_list = 1 : f (singleton z) [1, z] zs where
f s q6s p6s'@(p6:p6s)
| m < p6 = m : f (union (fromList $ map (* m) ps) s') q6s p6s'
| otherwise = f (union (fromList $ map (* p6) q6s) s) (p6:q6s) p6s
where ps = a027748_row m
(m, s') = deleteFindMin s
(z:zs) = a030516_list
-- Reinhard Zumkeller, Jun 03 2015
(Python)
from math import gcd
from sympy import factorint, integer_nthroot
def A069493(n):
def f(x):
c = n+x
for y1 in range(1, integer_nthroot(x, 11)[0]+1):
if all(d<=1 for d in factorint(y1).values()):
for y2 in range(1, integer_nthroot(z2:=x//y1**11, 10)[0]+1):
if gcd(y2, y1)==1 and all(d<=1 for d in factorint(y2).values()):
for y3 in range(1, integer_nthroot(z3:=z2//y2**10, 9)[0]+1):
if gcd(y3, y1)==1 and gcd(y3, y2)==1 and all(d<=1 for d in factorint(y3).values()):
for y4 in range(1, integer_nthroot(z4:=z3//y3**9, 8)[0]+1):
if gcd(y4, y1)==1 and gcd(y4, y2)==1 and gcd(y4, y3)==1 and all(d<=1 for d in factorint(y4).values()):
for y5 in range(1, integer_nthroot(z5:=z4//y4**8, 7)[0]+1):
if gcd(y5, y1)==1 and gcd(y5, y2)==1 and gcd(y5, y3)==1 and gcd(y5, y4)==1 and all(d<=1 for d in factorint(y5).values()):
c -= integer_nthroot(z5//y5**7, 6)[0]
return c
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
return bisection(f, n, n) # Chai Wah Wu, Sep 10 2024
CROSSREFS
KEYWORD
easy,nonn
AUTHOR
Benoit Cloitre, Apr 15 2002
STATUS
approved