login
a(1) = 1; for n>1, a(n) = smallest divisor d of a(n-1) not yet present in the sequence. If no such d exists, a(n) = 2*A007947(a(n-1)) + 1.
0

%I #10 Dec 04 2021 13:00:04

%S 1,3,7,15,5,11,23,47,95,19,39,13,27,9,7,15,31,63,21,43,87,29,59,119,

%T 17,35,71,143,287,41,83,167,335,67,135,45,31,63,43,87,175,25,11,23,47,

%U 95,191,383,767,1535,307,615,123,247,495,33,67,135,31,63,43

%N a(1) = 1; for n>1, a(n) = smallest divisor d of a(n-1) not yet present in the sequence. If no such d exists, a(n) = 2*A007947(a(n-1)) + 1.

%C Enters a loop starting at n = 420: [211, 423, 283, 567, 43, 87, 175, 71, 143, 287, 575, 231, 463, 927, 619, 1239, 2479, 4959, 3307, 6615].

%t m = d[1] = 1; {1}~Join~Reap[Do[Set[n, If[IntegerQ@ #1, #1, 1 + 2 SelectFirst[Reverse@ #2, SquareFreeQ]]] & @@ {SelectFirst[#, ! IntegerQ[d[#]] &], #} &@ Divisors[m]; Sow[n]; Set[d[n], i]; m = n, {i, 2, 61}]][[-1, -1]] (* _Michael De Vlieger_, Nov 09 2021 *)

%o (Python)

%o from sympy import primefactors, prod, divisors

%o terms = [1]

%o for i in range(100):

%o for j in divisors(terms[-1]):

%o if j not in terms:

%o terms.append(j)

%o break

%o else:

%o terms.append(prod(primefactors(terms[-1]))*2+1)

%o print(terms)

%Y Cf. A348871, A348872, A348881, A348873, A007947.

%K nonn

%O 1,2

%A _Gleb Ivanov_, Nov 09 2021