Lucas-Lehmer test: Difference between revisions

added ruby
(added Fortran)
(added ruby)
Line 7:
 
=={{header|Ada}}==
<ada>with Ada.Text_Io; use Ada.Text_Io;
with Ada.Integer_Text_Io; use Ada.Integer_Text_Io;
 
procedure Lucas_Lehmer_Test is
type Ull is mod 2**64;
function Mersenne(Item : Integer) return Boolean is
S : Ull := 4;
MP : Ull := 2**Item - 1;
begin
if Item = 2 then
return True;
else
for I in 3..Item loop
S := (S * S - 2) mod MP;
end loop;
return S = 0;
end if;
end Mersenne;
Upper_Bound : constant Integer := 64;
begin
Put_Line(" Mersenne primes:");
for P in 2..Upper_Bound loop
if Mersenne(P) then
Put(" M");
Put(Item => P, Width => 1);
end if;
end loop;
end Lucas_Lehmer_Test;</ada>
Output:
Mersenne primes:
Line 89:
{{works with|C99}}
Compiler options: gcc -std=c99 -lm Lucas-Lehmer_test.c -o Lucas-Lehmer_test<br>
<c>#include <math.h>
#include <stdio.h>
#include <limits.h>
#pragma precision=log10l(ULLONG_MAX)/2
 
typedef enum { FALSE=0, TRUE=1 } BOOL;
 
BOOL is_prime( int p ){
if( p == 2 ) return TRUE;
else if( p <= 1 || p % 2 == 0 ) return FALSE;
else {
BOOL prime = TRUE;
const int to = sqrt(p);
int i;
for(i = 3; i <= to; i+=2)
if (!(prime = p % i))break;
return prime;
}
}
 
BOOL is_mersenne_prime( int p ){
if( p == 2 ) return TRUE;
else {
const long long unsigned m_p = ( 1LLU << p ) - 1;
long long unsigned s = 4;
int i;
for (i = 3; i <= p; i++){
s = (s * s - 2) % m_p;
}
return s == 0;
}
}
 
int main(int argc, char **argv){
 
const int upb = log2l(ULLONG_MAX)/2;
int p;
printf(" Mersenne primes:\n");
for( p = 2; p <= upb; p += 1 ){
if( is_prime(p) && is_mersenne_prime(p) ){
printf (" M%u",p);
}
}
printf("\n");
}</java>
}
 
Output:
Line 142:
{{libheader|GMP}}
 
<precpp>#include <iostream>
#include <gmpxx.h>
 
Line 198:
//std::cout << "pow2( " << exp << " ) = " << ret << std::endl;
return ret;
}</precpp>
}
</pre>
 
Output: (Incomplete; It takes a long time.)
Line 255 ⟶ 254:
We use arbitrary-precision integers in order to be able to test any arbitrary prime.
 
<java>import java.math.BigInteger;
public class Mersenne
{
 
public static boolean isPrime(int p) {
if (p == 2)
return true;
else if (p <= 1 || p % 2 == 0)
return false;
else {
int to = (int)Math.sqrt(p);
for (int i = 3; i <= to; i += 2)
if (p % i == 0)
return false;
return true;
}
}
 
public static boolean isMersennePrime(int p) {
if (p == 2)
return true;
else {
BigInteger m_p = BigInteger.ONE.shiftLeft(p).subtract(BigInteger.ONE);
BigInteger s = BigInteger.valueOf(4);
for (int i = 3; i <= p; i++)
s = s.multiply(s).subtract(BigInteger.valueOf(2)).mod(m_p);
return s.equals(BigInteger.ZERO);
}
}
 
// an arbitrary upper bound can be given as an argument
public static void main(String[] args) {
int upb;
if (args.length == 0)
upb = 500;
else
upb = Integer.parseInt(args[0]);
 
System.out.println(" Finding Mersenne primes in M[2.." + upb + "]: ");
for (int p = 3; p <= upb; p += 2)
if (isPrime(p) && isMersennePrime(p))
System.out.print(" M" + p);
System.out.println();
}
}</java>
}
Output (after about eight hours):
Finding Mersenne primes in M[2..2147483647]:
Line 416 ⟶ 415:
 
=={{header|Python}}==
<python> from sys import stdout
from math import sqrt, log
 
def is_prime ( p ):
if p == 2: return True
elif p <= 1 or p % 2 == 0: return False
else:
for i in range(3, int(sqrt(p))+1, 2 ):
if p % i == 0: return False
return True
 
def is_mersenne_prime ( p ):
if p == 2:
return True
else:
m_p = ( 1 << p ) - 1;
s = 4;
for i in range(3, p+1):
s = (s ** 2 - 2) % m_p
return s == 0
 
precision = 20000 # maximum requested number of decimal places of 2 ** MP-1 #
long_bits_width = precision / log(2) * log(10)
upb_prime = int( long_bits_width - 1 ) / 2 # no unsigned #
upb_count = 45 # find 45 mprimes if int was given enough bits #
 
print (" Finding Mersenne primes in M[2..%d]:"%upb_prime)
 
count=0
for p in range(2, upb_prime+1):
if is_prime(p) and is_mersenne_prime(p):
print("M%d"%p),
stdout.flush()
count += 1
if count >= upb_count: break
print</python>
 
Output:
<pre> Finding Mersenne primes in M[2..33218]:
M2 M3 M5 M7 M13 M17 M19 M31 M61 M89 M107 M127 M521 M607 M1279 M2203 M2281 M3217 M4253 M4423 M9689 M9941 M11213 M19937 M21701 M23209</pre>
 
=={{header|Ruby}}==
<python>def is_prime ( p )
if p == 2
return true
elsif p <= 1 or p % 2 == 0
return false
else
(3 .. Math.sqrt(p)).step(2) do |i|
if p % i == 0
return false
end
end
return true
end
end
 
def is_mersenne_prime ( p )
if p == 2
return true
else
m_p = ( 1 << p ) - 1
s = 4
for i in 3..p
s = (s ** 2 - 2) % m_p
end
return s == 0
end
end
 
precision = 20000 # maximum requested number of decimal places of 2 ** MP-1 #
long_bits_width = precision / Math.log(2) * Math.log(10)
upb_prime = (long_bits_width - 1).to_i / 2 # no unsigned #
upb_count = 45 # find 45 mprimes if int was given enough bits #
 
puts " Finding Mersenne primes in M[2..%d]:"%upb_prime
 
count = 0
for p in 2..upb_prime
if is_prime(p) && is_mersenne_prime(p)
print "M%d "%p
count += 1
end
if count >= upb_count
break
end
end
puts</python>
 
Output:
Line 458 ⟶ 509:
 
=={{header|Scheme}}==
<scheme>;;;The heart of the algorithm
(define (S n)
(let ((m (- (expt 2 n) 1)))
(let loop ((c (- n 2)) (a 4))
(if (zero? c)
a
(loop (- c 1) (remainder (- (* a a) 2) m))))))
 
(define (mersenne-prime? n)
(if (= n 2)
#t
(zero? (S n))))
 
;;;Trivial unoptimized implementation for the base primes
(define (next-prime x)
(if (prime? (+ x 1))
(+ x 1)
(next-prime (+ x 1))))
 
(define (prime? x)
(let loop ((c 2))
(cond ((>= c x) #t)
((zero? (remainder x c)) #f)
(else (loop (+ c 1))))))
 
;;Main loop
(let loop ((i 45) (p 2))
(if (not (zero? i))
(if (mersenne-prime? p)
(begin
(display "M") (display p) (display " ")
(loop (- i 1) (next-prime p)))
(loop i (next-prime p)))))</scheme>
 
M2 M3 M5 M7 M13...
Anonymous user