# Algorithm Implementation/Mathematics/Primality Testing

In practice, primality testing for numbers of a size suitable for cryptographic applications has to be done in a probabilistic way. Such an algorithm can tell whether a given number is prime with extremely high probability, but cannot provide a certain proof. Many such algorithms operate with a variable number of rounds of testing; each additional round of testing makes it less likely to falsely classify the number, so the implementor can choose a suitable number of rounds to provide a trade off between the level of certainty and the running time of the algorithm.

## Pseudocode

```Input: n > 3, an odd integer to be tested for primality;
Input: k, a parameter that determines the accuracy of the test
Output: composite if n is composite, otherwise probably prime
write n − 1 as 2s·d with d odd by factoring powers of 2 from n − 1
LOOP: repeat k times:
pick a randomly in the range [2, n − 2]
if x = 1 or x = n − 1 then do next LOOP
for r = 1 .. s − 1
x ← x2 mod n
if x = 1 then return composite
if x = n − 1 then do next LOOP
return composite
return probably prime
```

## Python

```def is_probable_prime(n, k = 7):
"""use Rabin-Miller algorithm to return True (n is probably prime)
or False (n is definitely composite)"""
if n < 6:  # assuming n >= 0 in all cases... shortcut small cases here
return [False, False, True, True, False, True][n]
elif n & 1 == 0:  # should be faster than n % 2
return False
else:
s, d = 0, n - 1
while d & 1 == 0:
s, d = s + 1, d >> 1
for a in random.sample(xrange(2, min(n - 2, sys.maxint)), min(n - 4, k)):
x = pow(a, d, n)
if x != 1 and x + 1 != n:
for r in xrange(1, s):
x = pow(x, 2, n)
if x == 1:
return False  # composite for sure
elif x == n - 1:
a = 0  # so we know loop didn't continue to end
break  # could be strong liar, try another a
if a:
return False  # composite if we reached end of this loop
return True  # probably prime if reached end of outer loop
```

## Perl

```use Math::BigInt;
use POSIX qw/INT_MAX/;

sub is_probable_prime(\$\$){
my (\$n,\$k) = (Math::BigInt->new(\$_[0])->babs(),\$_[1]);
return 1 if \$n == 2 || \$n == 3;
return 0 if \$n < 5 || (\$n & 1) == 0;
my (\$s,\$d) = (0,\$n-1);
(\$s,\$d) = (\$s+1,\$d>>1) while (\$d & 1) == 0;
my \$maxa = (\$n-1 < int(INT_MAX) ? \$n-1 : int(INT_MAX));
for my \$i (0..\$k){
my \$a = Math::BigInt->new(int(rand(\$maxa-2)+2));
my \$x = \$a->bmodpow(\$d,\$n);
if(\$x != 1 && \$x+1 != \$n){
for my \$j (1..\$s){
\$x = \$x->bmodpow(2,\$n);
return 0 if \$x == 1;
if (\$x == \$n - 1){
\$a = 0;
last;
}
}
return 0 if \$a;
}
}
return 1;
}
```

## Java

```import java.math.BigInteger;
public class MillerRabin {

public static final BigInteger ZERO = BigInteger.ZERO;  // declaring constants
public static final BigInteger ONE  = BigInteger.ONE;
public static final BigInteger TWO = BigInteger.valueOf(2);
public static final int[] aValues = {2,3,5,7,11,13,17,19,23,29,31,37,41};  // values for bases

public static boolean testPr (BigInteger n, BigInteger a, int s, BigInteger d){
for (int i=0; i<s; i++){
BigInteger exp = TWO.pow(i);
exp = exp.multiply(d);
BigInteger res = a.modPow(exp, n);
if (res.equals(n.subtract(ONE)) || res.equals(ONE)){
return true;
}
}
return false;
}

public static boolean millerRabin (BigInteger n, int numValues){ // n = number to test
BigInteger d = n.subtract(ONE);                          // numValues = # of bases to test
int s = 0;
while (d.mod(TWO).equals(ZERO)){
s++;
d = d.divide(TWO);
}
System.out.print("Base ");
for (int i=0; i<numValues; i++){     // loops through the bases, terminating early if
BigInteger a = BigInteger.valueOf(aValues[i]);  // composite
boolean r = testPr(n, a, s, d);
System.out.print(aValues[i] + " ");
if (!r){
return false;
}
}
return true;
}
}
```