Yellowstone sequence

Yellowstone sequence is a draft programming task. It is not yet considered ready to be promoted as a complete task, for reasons that should be found in its talk page.

The Yellowstone sequence, also called the Yellowstone permutation, is defined as:

For n <= 3,

a(n) = n.

For n >= 4,

a(n) = the smallest number not already in sequence such that a(n) is relatively prime to a(n-1) and is not relatively prime to a(n-2).

The sequence is a permutation of the natural numbers, and gets its name from what its authors felt was a spiking, geyser like appearance of a plot of the sequence.

Task

Find and show as output the first 30 Yellowstone numbers.

Extra

Demonstrate how to plot, with x = n and y coordinate a(n), the first 100 Yellowstone numbers.

Example

a(4) is 4 because 4 is the smallest number following 1, 2, 3 in the sequence that is relatively prime to the entry before it (3), and is not relatively prime to the number two entries before it (2).

Related tasks

•   Greatest common divisor.
•   Plot coordinate pairs.

See also:

•   The OEIS entry:   A098550 The Yellowstone permutation.
•   Applegate et al, 2015: The Yellowstone Permutation [1]

Julia

<lang julia>using Plots

function yellowstone(N)

   a = [1, 2, 3]
   b = Dict(1 => 1, 2 => 1, 3 => 1)
   while length(a) < N
       for i in 4:typemax(Int)
           if !haskey(b, i) && (gcd(i, a[end]) == 1) && (gcd(i, a[end - 1]) > 1)
               push!(a, i)
               b[i] = 1
               break
           end
       end
   end
   return a

end

println("The first 30 entries of the Yellowstone permutation:\n", yellowstone(30))

x = 1:100 y = yellowstone(100) plot(x, y)

</lang>

Output:

The first 30 entries of the Yellowstone permutation:
[1, 2, 3, 4, 9, 8, 15, 14, 5, 6, 25, 12, 35, 16, 7, 10, 21, 20, 27, 22, 39, 11, 13, 33, 26, 45, 28, 51, 32, 17]

REXX

<lang rexx>/*REXX program calculates any number of terms in the Yellowstone (permutation) sequence.*/ parse arg m . /*obtain optional argument from the CL.*/ if m== | m=="," then m= 30 /*Not specified? Then use the default.*/ !.= 0 /*initialize an array of numbers(used).*/

1. = 0 /*count of Yellowstone numbers in seq. */

$= /*list " " " " " */

     do j=1  until #==m;  prev= # - 1
     if j<5  then do;  #= #+1;   @.#= j;  !.#= j;  !.j= 1;  $= strip($ j);  iterate;  end

        do k=1;   if !.k  then iterate          /*Already used?  Then skip this number.*/
        if gcd(k, @.#)\==1  |  gcd(k, @.prev)<2  then iterate   /*not meet requirement?*/
        #= #+1;   @.#= k;     !.k= 1;   $= $ k  /*bump ctr; assign; mark used; add list*/
        leave                                   /*find the next Yellowstone seq. number*/
        end   /*k*/
     end      /*j*/

say $ /*display a list of a Yellowstone seq. */ exit /*stick a fork in it, we're all done. */ /*──────────────────────────────────────────────────────────────────────────────────────*/ gcd: parse arg x,y; do until y==0; parse value x//y y with y x; end; return x</lang>

output   when using the default input:

1 2 3 4 9 8 15 14 5 6 25 12 35 16 7 10 21 20 27 22 39 11 13 33 26 45 28 51 32 17

zkl

Translation of: Julia

This sequence is limited to the max size of a Dictionary, 64k <lang zkl>fcn yellowstoneW{

  Walker.zero().tweak(fcn(a,b){
     foreach i in ([1..]){
        if(not b.holds(i) and i.gcd(a[-1])==1 and i.gcd(a[-2]) >1){

a.del(0).append(i); // only keep last two terms b[i]=True; return(i); }

     }
  }.fp(List(2,3), Dictionary(1,True, 2,True, 3,True))).push(1,2,3);

}</lang> <lang zkl>println("The first 30 entries of the Yellowstone permutation:"); yellowstoneW().walk(30).concat(", ").println();</lang>

Output:

The first 30 entries of the Yellowstone permutation:
1, 2, 3, 4, 9, 8, 15, 14, 5, 6, 25, 12, 35, 16, 7, 10, 21, 20, 27, 22, 39, 11, 13, 33, 26, 45, 28, 51, 32, 17

Plot using Gnuplot <lang zkl>gnuplot:=System.popen("gnuplot","w"); gnuplot.writeln("plot '-'"); yellowstoneW().pump(100,gnuplot.writeln,fcn(n){ String(" ",n) }); gnuplot.writeln("e"); gnuplot.flush(); ask("Hit return to finish"); gnuplot.close();</lang>