Problem Statement

Lucas-Lehmer Test: for

p

{\displaystyle p}

an odd prime, the Mersenne number

2

p

− 1

{\displaystyle 2^{p}-1}

is prime if and only if

2

p

− 1

{\displaystyle 2^{p}-1}

divides

S ( p − 1 )

{\displaystyle S(p-1)}

where

S ( n + 1 )

( S ( n )

)

2

− 2

{\displaystyle S(n+1)=(S(n))^{2}-2}

, and

S ( 1 )

4

{\displaystyle S(1)=4}

.

Calculate all Mersenne primes up to the implementation's maximum precision, or the 47th Mersenne prime   (whichever comes first).

Let's start with the solution:

  [ dup temp put
    dup bit 1 -
    4
    rot 2 - times
      [ dup *
        dup temp share >>
        dip [ over & ] +
        2dup > not if
          [ over - ]
        2 - ] 
    0 = 
    nip temp release ]    is l-l ( n --> b )

  25000 eratosthenes
  [] 25000 times [ i^ isprime if [ i^ join ] ]
  1 split
  witheach
    [ dup l-l iff join else drop ]
  echo