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}

.



Task

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




readlib 'mapm';
 
mapm.xdigits(100);
 
mersenne := proc(p::number) is
   local s, m;
   s := 4;
   m := mapm.xnumber(2^p) - 1;
   if p = 2 then
      return true
   else
      for i from 3 to p do
         s := (mapm.xnumber(s)^2 - 2) % m
      od;
      return mapm.xtoNumber(s) = 0
   fi
end;
 
for i from 3 to 64 do
   if mersenne(i) then
      write('M' & i & ' ')
   fi
od;