You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
 
 
 
 
 

82 lines
1.9 KiB

// Check whether a mersenne number is prime,
// using the Lucas-Lehmer test.
// [Donald Ervin Knuth: The Art of Computer Programming, Vol. II:
// Seminumerical Algorithms, second edition. Section 4.5.4, p. 391.]
// We work with integers.
#include <cl_integer.h>
// Checks whether 2^q-1 is prime, q an odd prime.
bool mersenne_prime_p (int q)
{
cl_I m = ((cl_I)1 << q) - 1;
int i;
cl_I L_i;
for (i = 0, L_i = 4; i < q-2; i++)
L_i = mod(L_i*L_i - 2, m);
return (L_i==0);
}
// Same thing, but optimized.
bool mersenne_prime_p_opt (int q)
{
cl_I m = ((cl_I)1 << q) - 1;
int i;
cl_I L_i;
for (i = 0, L_i = 4; i < q-2; i++) {
L_i = square(L_i) - 2;
L_i = ldb(L_i,cl_byte(q,q)) + ldb(L_i,cl_byte(q,0));
if (L_i >= m)
L_i = L_i - m;
}
return (L_i==0);
}
// Now we work with modular integers.
#include <cl_modinteger.h>
// Same thing, but using modular integers.
bool mersenne_prime_p_modint (int q)
{
cl_I m = ((cl_I)1 << q) - 1;
cl_modint_ring R = cl_find_modint_ring(m); // Z/mZ
int i;
cl_MI L_i;
for (i = 0, L_i = R->canonhom(4); i < q-2; i++)
L_i = R->minus(R->square(L_i),R->canonhom(2));
return R->equal(L_i,R->zero());
}
#include <cl_io.h> // we do I/O
#include <stdlib.h> // declares exit()
#include <cl_timing.h>
int main (int argc, char* argv[])
{
if (!(argc == 2)) {
fprint(cl_stderr, "Usage: lucaslehmer exponent\n");
exit(1);
}
int q = atoi(argv[1]);
if (!(q >= 2 && ((q % 2)==1))) {
fprint(cl_stderr, "Usage: lucaslehmer q with q odd prime\n");
exit(1);
}
bool isprime;
{ CL_TIMING; isprime = mersenne_prime_p(q); }
{ CL_TIMING; isprime = mersenne_prime_p_opt(q); }
{ CL_TIMING; isprime = mersenne_prime_p_modint(q); }
fprint(cl_stdout, "2^");
fprintdecimal(cl_stdout, q);
fprint(cl_stdout, "-1 is ");
if (isprime)
fprint(cl_stdout, "prime");
else
fprint(cl_stdout, "composite");
fprint(cl_stdout, "\n");
}
// Computing time on a i486, 33 MHz:
// 1279: 2.02 s
// 2281: 8.74 s
// 44497: 14957 s