Sylvester primes

Question

Sylvester's sequence can be defined recursively S(n) = S(n-1)*(S(n-1) + 1) for n >= 1 starting S(0) = 1.

Since S(n) and S(n) + 1 have no common divisors, it follows that S(n) has at least one more prime factor than S(n-1), and thus by induction, S(n) has at least n distinct prime factors.

This simple and constructive form of Euclid's proof of the infinity of primes was formulated by Filip Saidak.

To generate the sequence, select the smallest among the prime factors of S(n) that has not yet been selected. We call the infinite sequence constructed this way the Sylvester primes.

The first few steps generating the sequence are sketched by:

  n,  selected,       factors of S(n)            a(n)
  [1] {},             {2}                     ->   2,
  [2] {2},            {2, 3}                  ->   3,
  [3] {2, 3},         {2, 3, 7}               ->   7,
  [4] {2, 3, 7},      {2, 3, 7, 43}           ->  43,
  [5] {2, 3, 7, 43},  {2, 3, 7, 13, 43, 139}  ->  13.

Summary:

Input: A positive integer n.

Output: The sequence a(1), ..., a(n) of the first n Sylvester primes.

Optional: Comment if the program can generate the sequence for n = 10 in less than a quarter of an hour.

This is code-golf, so each language's shortest code in bytes wins.

Answer 1

JavaScript (ES11), 75 bytes

Expects a Bigint \$n\$ and returns \$a(n)\$ as another Bigint.

Computes \$a(1)\$ to \$a(12)\$ almost instantly.

f=(n,x=2n,h=x=>x%++v?h(x):h[v]?h(x/v--):--n)=>h(x,v=1n)?f(n,x*-~x,h[v]=h):v

Try it online!

(or 73 bytes if we assume that Sylvester's sequence is square-free, which is an unresolved problem)

Commented

f = (            // f is a recursive function taking:
  n,             //   n = input
  x = 2n,        //   x = current term of Sylvester's sequence
  h = x =>       //   h is a recursive helper function taking x
                 //   and setting v to the smallest prime factor
                 //   of x which was not listed before:
  x % ++v ?      //     increment v; if v is not a divisor of x:
    h(x)         //       try again with x unchanged
  :              //     else:
    h[v] ?       //       if h[v] is already defined:
      h(x / v--) //         try again with x/v and v-1
    :            //       else:
      --n        //         done: decrement n and return it
) =>             //
  h(x, v = 1n)   // invoke h, starting with v=1
  ?              // if n is not zero:
    f(           //   recursive call:
      n,         //     pass n
      x * -~x,   //     update x to the next term of the sequence
      h[v] = h   //     pass h with h[v] set
    )            //   end of recursive call
  :              // else:
    v            //   end of recursion: return v

---

Limited recursion, 82 bytes

This version can reach \$a(19)\$ within TIO's time limit.

f=(n,x=2n,h=y=>{for(v=1n;y%++v||h[y/=v,v]&&v--;);})=>h(x)||--n?f(n,x*-~x,h[v]=h):v

Try it online!

Answer 2

Jelly,  15  14 bytes

^{I'm not convinced this is as terse as possible.}

1‘×$СÆfḟṂṭɗ@/

A full program that accepts a positive integer, \$n\$, on stdin and prints a Jelly representation of the list of the first \$n\$ Sylvester primes.

Try it online!

How?

1‘×$СÆfḟṂṭɗ@/ - Main Link: no arguments
1              - literal one (first value of Current, below)
    С         - collect, repeating...
               - ...times: positive integer, N, from stdin (implicit)
   $           - ...action: last two links as a monad - f(Current):
 ‘             -              increment {Current} -> Current + 1
  ×            -              multiply {that] by {Current} -> Current × (Current + 1)
      Æf       - prime factors (vectorises)
             / - reduce by:
            @  -   with swapped arguments:
           ɗ   -     last three links as a dyad - f(right, left):
        ḟ      -       {right} filter discard {left}
         Ṃ     -       minimum of {that}
          ṭ    -       tack to {left}

Answer 3

05AB1E, 15 12 bytes

$FD>*Df¯Kß=ˆ

Given \$n\$, it'll output the first \$n\$ values on separated lines.

Try it online.

Original 15 bytes answer:

λD>*}fÅ»sDŠKß=ª

Outputs the infinite sequence, each value on a separated line.

Try it online.

Explanation:

$         # Push 1 and the input-integer
 F        # Pop and loop the input amount of times:
  D       #  Duplicate the current n
   >      #  Increase the copy by 1
    *     #  Multiply the two together: n*(n+1)
     D    #  Duplicate it for the next iteration
  f       #  Pop the copy, and push a list of its unique prime factors
   ¯      #  Push the global array of all previous terms
    K     #  Remove those values from the prime factors
     ß    #  Pop and leave the minimum
      =   #  Output it with trailing newline (without popping)
       ˆ  #  Pop and add it to the global array

λD>*}     # Generate the infinite Sylvester's sequence:
λ         #  Start the recursive environment
          #  to result in an infinite sequence
          #  Starting implicitly with a(0)=1
          #  And where every following a(n) is calculated by:
          #   (implicitly push the previous a(n-1))
 D        #   Duplicate it
  >       #   Increase the copy by 1
   *      #   Multiply the two together: a(n-1)*(a(n-1)+1)
    }     #  Close the recursive environment
f         # Map each inner value to a list of its unique prime factors
Å»        # Cumulative left-reduce this list by:
  s       #  Swap so the current list is at the top of the stack
   D      #  Duplicate it
    Š     #  Triple-swap it
          #  (the stack is now: list,primeFactors,list)
     K    #  Remove all values we've already output from the prime-factors list
      ß   #  Pop and leave the minimum remaining prime-factor
       =  #  Output it with trailing newline (without popping)
        ª #  Add it to the list for the next iteration

Answer 4

Stax, 15 bytes

^{I'm not convinced this is as terse as possible.}

Çö╫k1,`♠y&∞ìσ▀╕

Run and debug it at staxlang.xyz! Link is to unpacked version with a breakpoint to avoid freezing your browser.

Generates the sequence infinitely, printing to standard output.

Unpacked (17 bytes)

1Wc^*c:Fx-hQ]x+Xd
1                    push initial value 1
 W                   forever:
  c^*                  n*(n+1)
     c:F               distinct prime factors in increasing order
        x-             remove all already generated
          hQ           take head and print
            ]x+        add to already generated
               Xd      save in register x and discard from stack

This abuses x a little bit, saving a byte by not initializing it to an empty array. The first + adds array [2] and integer 0, and all the rest concatenate arrays. That means an extra 0 in the list, which is not a problem.

Answer 5

Charcoal, 27 bytes

≔¹ηＦＮ«≧×⊕ηη⊞υ²Ｗ﹪ηΠυ⊞υ⊕⊟υ»Ｉυ

Try it online! Link is to verbose version of code. Can calculate up to n=18 on TIO. Explanation: Assumes without proof that Sylvester's sequence is square-free.

≔¹η

Start with S(0)=1.

ＦＮ«

Repeat n times.

≧×⊕ηη

Calculate the next term of S.

⊞υ²

Start with 2 as the next candidate prime.

Ｗ﹪ηΠυ

Until the current S term is divisible by the product of the primes so far...

⊞υ⊕⊟υ

... increment the candidate.

»Ｉυ

Output the primes.

Answer 6

Ruby, 72 70 bytes

->n{r,*t=1;n.times{t<<(2..r*=r+1).find{|x|r%x<1&&t.all?{|y|x%y>0}}};t}

Try it online!

Calculates f[20] in 30 seconds on TIO.

Answer 7

Python, 185 173 bytes

-12 bytes thanks to @Sophia Antipolis

r=range
S=lambda j:j and S(j-1)*(S(j-1)+1)or 1
Q=lambda k:[x for x in r(1,S(k)+1)if S(k)%x<1and all(x%i for i in r(2,x//2))]
F=lambda i:[min({*Q(z+1)}-{*Q(z)})for z in r(i)]

Attempt This Online!

I think I'm a glutton for punishment, trying to do these with Python... It's fun, though!

• S is getting a list of Sylvester's sequence, recursively.
• Q is getting all of the primes of S
• F is getting a list of the minimum values for range 1 to i+1 (z) of Q(z), excluding values in Q(z-1)

This is certainly not a fast solution.

Answer 8

Maple, 105 bytes

For n = 7 returns {2, 3, 7, 13, 43, 139, 547} in 0.5 seconds.

a:=proc(n)p:=1;w:={};for k to n do p:=p^2+p;{seq(x[1],x=ifactors(p)[2])}minus w;w:=w union{min(%)}od end:

Answer 9

Perl 5 + Math::Prime::Util, 84 bytes

sub S{$_[0]?S($_[0]-1)*(S(-1+pop)+1):1}say+(grep!$s{$_}++,factor S($_))[0]for 1..pop

Try it online!

For input 10 it outputs 2, 3, 7, 43, 13, 3263443, 547, 29881, 5295435634831, 181. Runtime for n=10 varies a lot. From 12 seconds and up.

Answer 10

Japt, 19 bytes

Outputs the first n terms. Produces inaccurate results for n>6 as we get into scientific notation & floating point inaccuracies in the terms of the Sylvester Sequence.

È+²}h2ì)mk rÈpYkX Î

Try it