login
a(n) = A065642(A065642(A019565(n))) / A019565(n).
3

%I #20 May 09 2021 09:51:25

%S 1,4,9,3,25,4,5,3,49,4,7,3,7,4,5,3,121,4,9,3,11,4,5,3,11,4,7,3,7,4,5,

%T 3,169,4,9,3,13,4,5,3,13,4,7,3,7,4,5,3,13,4,9,3,11,4,5,3,11,4,7,3,7,4,

%U 5,3,289,4,9,3,17,4,5,3,17,4,7,3,7,4,5,3,17,4,9,3,11,4,5,3,11,4,7,3,7,4,5,3,17,4,9,3,13,4,5,3,13,4,7,3,7,4,5,3

%N a(n) = A065642(A065642(A019565(n))) / A019565(n).

%C 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.

%C 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)).

%H Antti Karttunen, <a href="/A285323/b285323.txt">Table of n, a(n) for n = 0..10000</a>

%F a(n) = A065642(A065642(A019565(n))) / A019565(n).

%o (PARI)

%o 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_

%o A007947(n) = factorback(factorint(n)[, 1]); \\ From _Andrew Lelechenko_, May 09 2014

%o A065642(n) = { my(r=A007947(n)); if(1==n,n,n = n+r; while(A007947(n) <> r, n = n+r); n); };

%o A285323(n) = A065642(A065642(A019565(n))) / A019565(n);

%o (Scheme)

%o (define (A285323 n) (/ (A065642 (A065642 (A019565 n))) (A019565 n)))

%o (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))))))

%o (Python)

%o from operator import mul

%o from sympy import prime, primefactors

%o from functools import reduce

%o 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_

%o def a007947(n): return 1 if n<2 else reduce(mul, primefactors(n))

%o def a065642(n):

%o if n==1: return 1

%o r=a007947(n)

%o n += r

%o while a007947(n)!=r:

%o n+=r

%o return n

%o def a(n): return a065642(a065642(a019565(n)))//a019565(n)

%o print([a(n) for n in range(101)]) # _Indranil Ghosh_, Apr 20 2017

%Y Cf. A000040, A007814, A014673, A020639, A019565, A065642, A285099, A285110, A285321, A285327.

%K nonn

%O 0,2

%A _Antti Karttunen_, Apr 19 2017