OFFSET
0,2
COMMENTS
After the initial a(0)=1, the third row of array A285321 divided by its first row. After 1, all terms are either primes or squares of primes. See A285110.
The sequence is completely determined by the positions of two least significant 1-bits of n: After initial zero, if n is a power of two (only one 1-bit present) or if prime(1+A285099(n)) > prime(1+A007814(n))^2, a(n) = prime(1+A007814(n))^2 = A020639(A019565(n))^2, otherwise a(n) = prime(1+A285099(n)) = A014673(A019565(n)).
LINKS
Antti Karttunen, Table of n, a(n) for n = 0..10000
PROG
(PARI)
A019565(n) = {my(j, v); factorback(Mat(vector(if(n, #n=vecextract(binary(n), "-1..1")), j, [prime(j), n[j]])~))}; \\ This function from M. F. Hasler
A007947(n) = factorback(factorint(n)[, 1]); \\ From Andrew Lelechenko, May 09 2014
(Scheme)
(define (A285323 n) (cond ((zero? n) 1) ((or (= 1 (A000120 n)) (> (A000040 (+ 1 (A285099 n))) (A000290 (A000040 (+ 1 (A007814 n)))))) (A000290 (A000040 (+ 1 (A007814 n))))) (else (A000040 (+ 1 (A285099 n))))))
(Python)
from operator import mul
from sympy import prime, primefactors
from functools import reduce
def a019565(n): return reduce(mul, (prime(i+1) for i, v in enumerate(bin(n)[:1:-1]) if v == '1')) if n > 0 else 1 # This function from Chai Wah Wu
def a007947(n): return 1 if n<2 else reduce(mul, primefactors(n))
def a065642(n):
if n==1: return 1
r=a007947(n)
n += r
while a007947(n)!=r:
n+=r
return n
def a(n): return a065642(a065642(a019565(n)))//a019565(n)
print([a(n) for n in range(101)]) # Indranil Ghosh, Apr 20 2017
CROSSREFS
KEYWORD
nonn
AUTHOR
Antti Karttunen, Apr 19 2017
STATUS
approved