Mersenne primes

Mersenne primes 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.

Mersenne primes:

Challenge:

Create code that will list (preferably calculate) all of the Mersenne primes until some limitation is reached. For information on what a Mersenne prime is, go to this link: [[1]]

AppleScript

<lang> on isPrime(integ) set isComposite to "" if (integ / 2) = (integ / 2 div 1) then log integ & " is composite because 2 is a factor" as string --buttons {"OK", "Cancel"} default button 1 cancel button 2

else set x to 2 set sqrtOfInteg to integ ^ 0.5 repeat until x = integ ^ 0.5 + 1 as integer if (integ / x) = integ / x div 1 then log integ & " is composite because " & x & " & " & (integ / x div 1) & " are factors" as string --buttons {"OK", "Cancel"} default button 1 cancel button 2 set isComposite to 1 set x to x + 1 else

set x to x + 1 end if

end repeat log integ & " is prime" as string --buttons {"OK", "Cancel"} default button 1 cancel button 2 if isComposite = 1 then log integ & "is composite" else display dialog integ end if end if

end isPrime set x to 2 repeat isPrime(((2 ^ x) - 1) div 1) set x to x + 1 end repeat </lang>

D

Simplest thing that could possibly work. Using better primality tests will allow for more results to be calculated in a reasonable amount of time. <lang D>import std.bigint; import std.stdio;

bool isPrime(BigInt bi) {

   if (bi < 2) return false;
   if (bi % 2 == 0) return bi == 2;
   if (bi % 3 == 0) return bi == 3;
   
   auto test = BigInt(5);
   while (test * test < bi) {
       if (bi % test == 0) return false;
       test += 2;
       if (bi % test == 0) return false;
       test += 4;
   }

   return true;

}

void main() {

   auto base = BigInt(2);

   for (int pow=1; pow<32; pow++) {
       if (isPrime(base-1)) {
           writeln("2 ^ ", pow, " - 1");
       }
       base *= 2;
   }

}</lang>

Output:

2 ^ 2 - 1
2 ^ 3 - 1
2 ^ 5 - 1
2 ^ 7 - 1
2 ^ 13 - 1
2 ^ 17 - 1
2 ^ 19 - 1
2 ^ 31 - 1

Go

The github.com/ncw/gmp package is a drop-in replacement for Go's math/big package. It's a CGo wrapper around the C GMP library and under these circumstances is two to four times as fast as the native Go package. Editing just the import line you can use whichever is more convenient for you (CGo has drawbacks, including limited portability). Normally build tags would be used to control this instead of editing imports in the source, but this keeps the example simpler. <lang go>package main

import ( "fmt" "time"

// Use one or the other of these: "math/big" //big "github.com/ncw/gmp" )

func main() { start := time.Now() one := big.NewInt(1) mp := big.NewInt(0) bp := big.NewInt(0) const max = 22 for count, p := 0, uint(2); count < max; { mp.Lsh(one, p) mp.Sub(mp, one) if mp.ProbablyPrime(0) { elapsed := time.Since(start).Seconds() if elapsed >= 0.01 { fmt.Printf("2 ^ %-4d - 1 took %6.2f secs\n", p, elapsed) } else { fmt.Printf("2 ^ %-4d - 1\n", p) } count++ } for { if p > 2 { p += 2 } else { p = 3 } bp.SetUint64(uint64(p)) if bp.ProbablyPrime(0) { break } } } }</lang>

Output:

Results using the GMP package on a 3.4 GHz Xeon E3-1245:

2 ^ 2    - 1
2 ^ 3    - 1
2 ^ 5    - 1
2 ^ 7    - 1
2 ^ 13   - 1
2 ^ 17   - 1
2 ^ 19   - 1
2 ^ 31   - 1
2 ^ 61   - 1
2 ^ 89   - 1
2 ^ 107  - 1
2 ^ 127  - 1
2 ^ 521  - 1
2 ^ 607  - 1
2 ^ 1279 - 1 took   0.05 secs
2 ^ 2203 - 1 took   0.38 secs
2 ^ 2281 - 1 took   0.44 secs
2 ^ 3217 - 1 took   1.53 secs
2 ^ 4253 - 1 took   4.39 secs
2 ^ 4423 - 1 took   5.02 secs
2 ^ 9689 - 1 took  73.78 secs
2 ^ 9941 - 1 took  81.24 secs

(A previous run on more modest hardware was ~365 seconds for M₉₉₄₁.)

Java

Translation of: Kotlin

<lang Java>import java.math.BigInteger;

public class MersennePrimes {

   private static final int MAX = 20;

   private static final BigInteger ONE = BigInteger.ONE;
   private static final BigInteger TWO = BigInteger.valueOf(2);

   private static boolean isPrime(int n) {
       if (n < 2) return false;
       if (n % 2 == 0) return n == 2;
       if (n % 3 == 0) return n == 3;
       int d = 5;
       while (d * d <= n) {
           if (n % d == 0) return false;
           d += 2;
           if (n % d == 0) return false;
           d += 4;
       }
       return true;
   }

   public static void main(String[] args) {
       int count = 0;
       int p = 2;
       while (true) {
           BigInteger m = TWO.shiftLeft(p - 1).subtract(ONE);
           if (m.isProbablePrime(10)) {
               System.out.printf("2 ^ %d - 1\n", p);
               if (++count == MAX) break;
           }
           // obtain next prime, p
           do {
               p = (p > 2) ? p + 2 : 3;
           } while (!isPrime(p));
       }
   }

}</lang>

Output:

2 ^ 2 - 1
2 ^ 3 - 1
2 ^ 5 - 1
2 ^ 7 - 1
2 ^ 13 - 1
2 ^ 17 - 1
2 ^ 19 - 1
2 ^ 31 - 1
2 ^ 61 - 1
2 ^ 89 - 1
2 ^ 107 - 1
2 ^ 127 - 1
2 ^ 521 - 1
2 ^ 607 - 1
2 ^ 1279 - 1
2 ^ 2203 - 1
2 ^ 2281 - 1
2 ^ 3217 - 1
2 ^ 4253 - 1
2 ^ 4423 - 1

Julia

Works with: Julia version 0.6

Julia module Primes uses Miller-Rabin primality test.

<lang julia>using Primes

mersenne(n::Integer) = convert(typeof(n), 2) ^ n - one(n) function main(nmax::Integer)

   n = ith = zero(nmax)
   while ith ≤ nmax
       if isprime(mersenne(n))
           println("M$n")
           ith += 1
       end
       n += 1
   end

end

main(big(20))</lang>

Output:

M2
M3
M5
M7
M13
M17
M19
M31
M61
M89
M107
M127
M521
M607
M1279
M2203
M2281
M3217
M4253
M4423
M9689

Kotlin

This task is similar to the Lucas-Lehmer test task except that you can use whatever method you like to test the primality of the Mersenne numbers. Here, I've chosen to use the JDK's BigInteger.isProbablePrime(certainty) method. The exact algorithm is implementation dependent --- GNU classpath uses only Miller-Rabin, while Oracle JDK uses Miller-Rabin and sometimes adds a Lucas test (this is not the Lucas-Lehmer test).

A 'certainty' parameter of 10 is enough to find the first 20 Mersenne primes but as even this takes about 90 seconds on my modest machine I've not bothered going beyond that. <lang scala>// version 1.2.10

import java.math.BigInteger

const val MAX = 20

val bigOne = BigInteger.ONE val bigTwo = 2.toBigInteger()

/* for checking 'small' primes */ fun isPrime(n: Int): Boolean {

   if (n < 2) return false
   if (n % 2 == 0) return n == 2
   if (n % 3 == 0) return n == 3
   var d : Int = 5
   while (d * d <= n) {
       if (n % d == 0) return false
       d += 2
       if (n % d == 0) return false
       d += 4
   }
   return true

}

fun main(args: Array<String>) {

   var count = 0
   var p = 2
   while (true) {
       val m = (bigTwo shl (p - 1)) - bigOne
       if (m.isProbablePrime(10)) {
           println("2 ^ $p - 1")
           if (++count == MAX) break
       }
       // obtain next prime, p
       while(true) {
           p = if (p > 2) p + 2 else 3
           if (isPrime(p)) break
       }
   }

}</lang>

Output:

2 ^ 2 - 1
2 ^ 3 - 1
2 ^ 5 - 1
2 ^ 7 - 1
2 ^ 13 - 1
2 ^ 17 - 1
2 ^ 19 - 1
2 ^ 31 - 1
2 ^ 61 - 1
2 ^ 89 - 1
2 ^ 107 - 1
2 ^ 127 - 1
2 ^ 521 - 1
2 ^ 607 - 1
2 ^ 1279 - 1
2 ^ 2203 - 1
2 ^ 2281 - 1
2 ^ 3217 - 1
2 ^ 4253 - 1
2 ^ 4423 - 1

PARI/GP

<lang parigp>LL(p)={

 my(m=Mod(4,1<<p-1));
 for(i=3,p,m=m^2-2);
 m==0

}; forprime(p=2,, if(LL(p), print("2^"p"-1")))</lang>

Perl

Since GIMPS went to the trouble of dedicating thousands of CPU years to finding Mersenne primes, we should be kind enough to use the results. The ntheory module front end does this, so the results up to 43 million is extremely fast (4 seconds), and we can reduce this another 10x by only checking primes. After the GIMPS double-checked mark, a Lucas-Lehmer test is done using code similar to Rosetta Code Lucas-Lehmer in C+GMP.

If this is too contrived, we can use Math::Prime::Util::GMP::is_mersenne_prime instead, which will run the Lucas-Lehmer test on each input. The first 23 Mersenne primes are found in under 15 seconds.

Library: ntheory

<lang perl>use ntheory qw/forprimes is_mersenne_prime/; forprimes { is_mersenne_prime($_) && say } 1e9;</lang>

Output:

2
3
5
7
13
17
19
31
61
...

Perl 6

Works with: Rakudo version 2018.01

We already have a multitude of tasks that demonstrate how to find Mersenne primes; Prime decomposition, Primality by trial division, Trial factoring of a Mersenne number, Lucas-Lehmer test, Miller–Rabin primality_test, etc. that all have Perl 6 entries. I'm not sure what I could add here that would be useful.

Hmmm.

Create code that will list all of the Mersenne primes until some limitation is reached.

It doesn't specify to calculate them, only to list them; why throw away all of the computer millenia of processing power that the GIMPS has invested?

<lang perl6>use HTTP::UserAgent; use Gumbo;

my $table = parse-html(HTTP::UserAgent.new.get('https://www.mersenne.org/primes/').content, :TAG

); say 'All known Mersenne primes as of ', Date(now); say 'M', ++$, ": 2$_ - 1" for $table[1]».[*][0][*].comb(/'exp_lo='\d+/)».subst(/\D/, ,:g) .trans([<0123456789>.comb] => [<⁰¹²³⁴⁵⁶⁷⁸⁹>.comb]).words; </lang>

Output:

All known Mersenne primes as of 2018-01-27
M1: 2² - 1
M2: 2³ - 1
M3: 2⁵ - 1
M4: 2⁷ - 1
M5: 2¹³ - 1
M6: 2¹⁷ - 1
M7: 2¹⁹ - 1
M8: 2³¹ - 1
M9: 2⁶¹ - 1
M10: 2⁸⁹ - 1
M11: 2¹⁰⁷ - 1
M12: 2¹²⁷ - 1
M13: 2⁵²¹ - 1
M14: 2⁶⁰⁷ - 1
M15: 2¹²⁷⁹ - 1
M16: 2²²⁰³ - 1
M17: 2²²⁸¹ - 1
M18: 2³²¹⁷ - 1
M19: 2⁴²⁵³ - 1
M20: 2⁴⁴²³ - 1
M21: 2⁹⁶⁸⁹ - 1
M22: 2⁹⁹⁴¹ - 1
M23: 2¹¹²¹³ - 1
M24: 2¹⁹⁹³⁷ - 1
M25: 2²¹⁷⁰¹ - 1
M26: 2²³²⁰⁹ - 1
M27: 2⁴⁴⁴⁹⁷ - 1
M28: 2⁸⁶²⁴³ - 1
M29: 2¹¹⁰⁵⁰³ - 1
M30: 2¹³²⁰⁴⁹ - 1
M31: 2²¹⁶⁰⁹¹ - 1
M32: 2⁷⁵⁶⁸³⁹ - 1
M33: 2⁸⁵⁹⁴³³ - 1
M34: 2¹²⁵⁷⁷⁸⁷ - 1
M35: 2¹³⁹⁸²⁶⁹ - 1
M36: 2²⁹⁷⁶²²¹ - 1
M37: 2³⁰²¹³⁷⁷ - 1
M38: 2⁶⁹⁷²⁵⁹³ - 1
M39: 2¹³⁴⁶⁶⁹¹⁷ - 1
M40: 2²⁰⁹⁹⁶⁰¹¹ - 1
M41: 2²⁴⁰³⁶⁵⁸³ - 1
M42: 2²⁵⁹⁶⁴⁹⁵¹ - 1
M43: 2³⁰⁴⁰²⁴⁵⁷ - 1
M44: 2³²⁵⁸²⁶⁵⁷ - 1
M45: 2³⁷¹⁵⁶⁶⁶⁷ - 1
M46: 2⁴²⁶⁴³⁸⁰¹ - 1
M47: 2⁴³¹¹²⁶⁰⁹ - 1
M48: 2⁵⁷⁸⁸⁵¹⁶¹ - 1
M49: 2⁷⁴²⁰⁷²⁸¹ - 1
M50: 2⁷⁷²³²⁹¹⁷ - 1

Ring

<lang ring>

1. Project : Mersenne primes
2. Date  : 2018/01/22
3. Author : Gal Zsolt [~ CalmoSoft ~]
4. Email  : <calmosoft@gmail.com>

n = 0 while true

       n = n +1
       if isprime(pow(2,n)-1) = 1
          see n + nl
       ok

end

func isprime num

      if (num <= 1) return 0 ok
      if (num % 2 = 0) and num != 2 return 0 ok
      for i = 3 to floor(num / 2) -1 step 2
           if (num % i = 0) return 0 ok
      next
      return 1

</lang> Output:

2
3
5
7
13
17
19

Scala

<lang Scala>object MersennePrimes extends App {

 import Stream._

 def primeSieve(s: Stream[Int]): Stream[Int] =
   s.head #:: primeSieve(s.tail filter { _ % s.head != 0 })
 val primes = primeSieve(from(2))

 def mersenne(p: Int): BigInt = (BigInt(2) pow p) - 1

 def s(mp: BigInt, p: Int): BigInt = { if (p == 1) 4 else ((s(mp, p - 1) pow 2) - 2) % mp }

 val upbPrime = 9941
 println(s"Finding Mersenne primes in M[2..$upbPrime]")
 ((primes takeWhile (_ <= upbPrime)).par map { p => (p, mersenne(p)) }
   map { p => if (p._1 == 2) (p, 0) else (p, s(p._2, p._1 - 1)) } filter { _._2 == 0 })
   .foreach { p =>
     println(s"prime M${(p._1)._1}: " +
       { if ((p._1)._1 < 200) (p._1)._2 else s"(${(p._1)._2.toString.size} digits)" })
   }
 println("That's All Folks!")

}</lang>

Sidef

Uses the is_mersenne_prime() function from Math::Prime::Util::GMP. <lang ruby>for p in (^Inf -> lazy.grep { .is_mersenne_prime }) {

   say "2^#{p} - 1"

}</lang>

Output:

2^2 - 1
2^3 - 1
2^5 - 1
2^7 - 1
2^13 - 1
2^17 - 1
2^19 - 1
2^31 - 1
2^61 - 1
2^89 - 1
2^107 - 1
2^127 - 1
2^521 - 1
2^607 - 1
2^1279 - 1
2^2203 - 1
2^2281 - 1
2^3217 - 1
2^4253 - 1
2^4423 - 1
2^9689 - 1
2^9941 - 1
^C
sidef mersenne.sf  12.47s user 0.02s system 99% cpu 12.495 total

zkl

Uses libGMP (GNU MP Bignum Library) and its Miller-Rabin probabilistic primality testing algorithm. <lang zkl>var [const] BN=Import("zklBigNum"); // libGMP fcn mprimes{

  n,m := BN(2),0;
  foreach e in ([2..]){
     n,m = n.shiftLeft(1), n-1;
     if(m.probablyPrime()) println("2^%d - 1".fmt(e));
  }

}() // gets rather slow after M(4423)</lang>

Output:

2^2 - 1
2^3 - 1
2^5 - 1
2^7 - 1
2^13 - 1
2^17 - 1
2^19 - 1
2^31 - 1
2^61 - 1
2^89 - 1
2^107 - 1
2^127 - 1
2^521 - 1
2^607 - 1
2^1279 - 1
2^2203 - 1
2^2281 - 1
2^3217 - 1
2^4253 - 1
2^4423 - 1
2^9689 - 1
2^9941 - 1
2^11213 - 1
...