Miller–Rabin primality test: Difference between revisions

Content added Content deleted
(Added a much more readable, simpler and faster algorithm.)
(Updated D entry)
Line 840: Line 840:
<lang d>import std.random;
<lang d>import std.random;


bool isProbablePrime(in ulong n, in int k) {
bool isProbablePrime(in ulong n, in uint k=10) /*nothrow*/ @safe /*@nogc*/ {
static long modPow(long b, long e, in long m)
static ulong modPow(ulong b, ulong e, in ulong m)
pure nothrow {
pure nothrow @safe @nogc {
long result = 1;
ulong result = 1;
while (e > 0) {
while (e > 0) {
if ((e & 1) == 1) {
if ((e & 1) == 1)
result = (result * b) % m;
result = (result * b) % m;
}
b = (b ^^ 2) % m;
b = (b * b) % m;
e >>= 1;
e >>= 1;
}
}
Line 854: Line 853:
}
}


if (n < 2 || n % 2 == 0)
if (n < 2 || n % 2 == 0)
return n == 2;
return n == 2;


Line 863: Line 862:
s++;
s++;
}
}
assert(2 ^^ s * d == n - 1);
assert(2 ^^ s * d == n - 1);


outer:
outer:
foreach (_; 0 .. k) {
foreach (immutable _; 0 .. k) {
ulong a = uniform(2, n);
immutable ulong a = uniform(2, n);
ulong x = modPow(a, d, n);
ulong x = modPow(a, d, n);
if (x == 1 || x == n - 1)
if (x == 1 || x == n - 1)
continue;
continue;
foreach (__; 1 .. s) {
foreach (immutable __; 1 .. s) {
x = modPow(x, 2, n);
x = modPow(x, 2, n);
if (x == 1) return false;
if (x == 1)
if (x == n - 1) continue outer;
return false;
if (x == n - 1)
continue outer;
}
}
return false;
return false;
Line 882: Line 883:
}
}


void main() { // demo code
void main() { // Demo code.
import std.stdio, std.range, std.algorithm;
import std.stdio, std.range, std.algorithm;

writeln(filter!(n => isProbablePrime(n, 10))(iota(2, 30)));
iota(2, 30).filter!isProbablePrime.writeln;
}</lang>
}</lang>
{{out}}
{{out}}