# Miller Rabin Probabilistic Primality Algorithm

Miller Rabin is undoubtedly the most popular probabilistic primality testing algorithm.

It uses the concept of Fermat’s Little Theorem.So if you don’t know fermat’s little theorem,it is advised to understand it first – Fermat

Let us Begin.

For an odd prime number p, we can write $\mathbf{\left ( p-1 \right )= {2^{d}.s}}$

Here s is odd and d>=1

By Fermat’s Little Theorem,we can write –

$\mathbf{a^{\left&space;(&space;p-1&space;\right&space;)}-1&space;\equiv&space;0&space;\pmod&space;p}$

$\mathbf{a^{2^{d}.s}-1&space;\equiv&space;0&space;\pmod&space;p}$

$\Rightarrow&space;\left&space;(&space;{a^{2^{d-1}.s}}&space;\right&space;)^2-1\equiv&space;0\pmod&space;p$

$\Rightarrow&space;\left&space;(&space;{a^{2^{d-1}.s}}&space;-1&space;\right&space;).\left&space;(&space;{a^{2^{d-1}.s}}&space;+1&space;\right&space;)\equiv&space;0\pmod&space;p$

$\Rightarrow&space;\left&space;(&space;{a^{2^{d-2}.s}}&space;-1&space;\right&space;).\left&space;(&space;{a^{2^{d-2}.s}}&space;+1&space;\right&space;).\left&space;(&space;{a^{2^{d-1}.s}}&space;+1&space;\right&space;)\equiv&space;0\pmod&space;p$

$\Rightarrow&space;\left&space;(&space;{a^{s}}&space;-1&space;\right&space;).\left&space;(&space;{a^{s}}&space;+1&space;\right&space;).\left&space;(&space;{a^{2s}}&space;+1&space;\right&space;).\left&space;(&space;{a^{4s}}&space;+1&space;\right&space;).......\left&space;(&space;{a^{2^{d-1}.s}}&space;+1&space;\right&space;)\equiv&space;0\pmod&space;p$

Thus for p to be prime , at least one of the terms in the expansion should be divisible by p.

In other words,

If none of the conditions are true then is definitely composite, otherwise p maybe prime.

To increase the accuracy of the algorithm, we will perform the check for certain number of iterations for some random values of a lesser than p.

Here’s  the python implementation –

import random

def exp(a,p,m):
if (p==0):
return 1
elif (p%2==1):
return (a*exp(a,p-1,m))%m
x=exp(a,p/2,m)
return (x*x)%m

def miller_rabin(p):

if (p==2 or p==3):return 1
if (p%2==0):return 0
if (p<3):return 0

for i in range(10):
#random 'a'
a = random.randrange(2,p-2)

s=p-1
while(s%2==0):s/=2

mod = exp(a,s,p)
if (mod==1 or mod==p-1):
continue

flag=0
while(s!=(p-1)):
mod = exp(a,s,p)
if (mod==(p-1)):
flag=1
break
s*=2
if (flag==0):
return 0
return 1


The Algorithm is slower than Fermat’s Little primality test but is still preferred due to it’s high accuracy.

NOTE – If you are implementing the algorithm in C/C++ , take care of the overflows while multiplying two long long ints  in modular exponentiation.

Try out this problem on SPOJ/PON

For any queries,drop a comment below 🙂