The OEIS mourns the passing of Jim Simons and is grateful to the Simons Foundation for its support of research in many branches of science, including the OEIS.
login
The OEIS is supported by the many generous donors to the OEIS Foundation.

 

Logo
Hints
(Greetings from The On-Line Encyclopedia of Integer Sequences!)
A358093 a(n) = n for 1 <= n <= 2; thereafter a(n) is the least unused m such that rad(m) = rad(rad(a(n-1)) + rad(a(n-2))), where rad(m) = A007947(m). 2

%I #65 Nov 14 2022 00:34:30

%S 1,2,3,5,4,7,9,10,13,23,6,29,35,8,37,39,38,77,115,12,11,17,14,31,15,

%T 46,61,107,42,149,191,170,19,21,20,961,41,18,47,53,40,63,29791,26,57,

%U 83,70,51,121,62,73,45,22,1369,59,24,65,71,34,105,139,122,87,209,74

%N a(n) = n for 1 <= n <= 2; thereafter a(n) is the least unused m such that rad(m) = rad(rad(a(n-1)) + rad(a(n-2))), where rad(m) = A007947(m).

%C In other words, a(n) is the least m having the same squarefree kernel as the sum of the squarefree kernels of a(n-1) and a(n-2).

%C Primes do not appear in natural order.

%C Conjectured to be a permutation of the positive integers.

%C From _Michael De Vlieger_, Nov 09 2022: (Start)

%C Let i = rad(a(n-2)), j = rad(a(n-1)), and s = rad(i+j). If s = p (prime), then a(n) = m = p^e, e >= 1. Because a(n) = m such that rad(m) = p, the j-th occasion of s = p implies a(n) = p^j.

%C By the same token, generally, the j-th occasion of s implies a(n) = M*s, where M is regular to s, i.e., M is a product restricted to primes p|s, with M appearing in order. For example, if s = 10, then the j-th occasion of s = 10 implies a(n) = 10*A003592(j). These are consequences of conservation of squarefree kernel s, the lexical and greedy nature of the sequence.

%C It is clear that a(n) >= s, since rad(a(n)) = s, hence s | a(n), and on account of the lexical and greedy nature of the sequence.

%C Prime a(n) implies s = a(n), while s < a(n) for a(n) that are composite prime powers (i.e., a(n) in A246547) since s is in those cases prime.(End)

%C From _Michael De Vlieger_, Nov 13 2022: (Start)

%C Adjacent terms i and j are coprime as consequence of s mod p = i mod p + j mod p and a(n) = M*s such that rad(M) = rad(s). p | a(n-2) and p | a(n-1) implies p | a(n), and subsequent adjacent terms are likewise divisible by p, contradicting initial a(1) coprime to a(2). Because the sequence starts with dissimilar residues (mod p), p divides only one of {a(n-2), a(n-1), a(n)}.

%C 2 | a(n) for n mod 3 = 2. (End)

%H Michael De Vlieger, <a href="/A358093/b358093.txt">Table of n, a(n) for n = 1..1500</a>

%H Michael De Vlieger, <a href="/A358093/a358093.png">Log log scatterplot of a(n)</a>, n = 1..1500, highlighting primes in red, other prime powers in gold, and other numbers in blue.

%e a(5) = 4 because rad(rad(3) + rad(5)) = rad(3 + 5) = rad(8) = 2, and 4 is the least unused number m such that rad(m) = 2.

%e a(36) = rad(rad(a(34)) + rad(a(35)) = rad(rad(21) + rad(20)) = rad(31) = 31, and since 31 has appeared at a(24), a(36) = 31^2 = 961.

%e a(43) = 29791 (31^3) because rad(40) + rad(63) = 31; we already have 31 and 31^2.

%t nn = 65; c[_] = False; p[_] = q[_] = 1; Array[Set[{a[#], c[#]}, {#, True}] &, 3]; Array[(q[#]++; p[#]++) &[Times @@ FactorInteger[a[#]][[All, 1]]] &, 2, 2]; u = 4; Set[{i, j}, Map[Times @@ FactorInteger[#][[All, 1]] &, {a[2], a[3]}]]; Do[k = Times @@ FactorInteger[i + j][[All, 1]]; If[PrimeQ[k], m = k^p[k]; p[k]++, m = q[k]; While[Nand[! c[m k], PowerMod[k, k, m] == 0], m++]; m *= k]; q[k]++; Set[{a[n], c[m], i, j}, {m, True, j, k}]; If[m == u, While[c[u], u++]], {n, 4, nn}]; Array[a, nn] (* _Michael De Vlieger_, Nov 09 2022 *)

%o (PARI) first(n) = { n = max(n, 3); res = vector(n); for(i = 1, 3, res[i] = i); for(i = 4, n, target = rad(rad(res[i-1]) + rad(res[i-2])); resset = Set(res); pr = factor(target)[, 1]; step = target; if(omega(target) >= 1, forstep(j = target, oo, step, if(rad(j) == target, if(#setminus(Set(j), resset) == 1, res[i] = j; next(2) ) ) ) , for(j = 1, oo, if(#setminus(Set(target ^ j), resset) == 1, res[i] = j; next(2) ) ) ) ); res }

%o rad(n) = factorback(factorint(n)[, 1]); \\ _David A. Corneth_, Nov 09 2022

%o (PARI) nn=300; a=List([1,2]); rad(n)= factorback(factorint(n)[, 1]);

%o {for(i=3,nn, listput(a,1); j=rad(rad(a[i-1])+rad(a[i-2])); while((rad(j)!=rad(rad(a[i-1])+rad(a[i-2])) || #select(n->n==j, a) > 0),j+= rad(rad(a[i-1])+rad(a[i-2]))); a[i]=j; ); print(a);} \\ _Gerry Martens_, Nov 10 2022

%Y Cf. A007947, A005117, A121369, A354184.

%K nonn

%O 1,2

%A _David James Sycamore_, Nov 08 2022

%E More terms from _David A. Corneth_, Nov 09 2022

Lookup | Welcome | Wiki | Register | Music | Plot 2 | Demos | Index | Browse | More | WebCam
Contribute new seq. or comment | Format | Style Sheet | Transforms | Superseeker | Recents
The OEIS Community | Maintained by The OEIS Foundation Inc.

License Agreements, Terms of Use, Privacy Policy. .

Last modified June 9 15:27 EDT 2024. Contains 373244 sequences. (Running on oeis4.)