#include <iostream>
#include <ctime>
using namespace std;
int RandInt(int Low, int High)        //routine to return a random number between
{                    //low and high boundaries
return rand( ) % ( High - Low + 1 ) + Low;
}

int Witness(int A, int i, int N)         //if witness does not return 1 then N is
{                                            //definately composite. Fermat’s little
int X, Y;                                 //theorem, says that if n is prime, then
//for all a between 1 and n - 1, an - 1 mod n = 1.
if( i == 0 )                              //the contrapositive is if an - 1 mod n ? 1,
return 1;                              //then n is not prime. Take a randomly!

X = Witness( A, i / 2, N );
if( X == 0 )                              // if N is recursively composite, stop
return 0;

Y = ( X * X ) % N;                    // N is not prime if we find a non-trivial root of 1
if( Y == 1 && X != 1 && X != N - 1 )
return 0;

if( i % 2 != 0 )
Y = ( A * Y ) % N;
return Y;
}

int CheckPrime( int N )                 //repeat this procedure as many times as
{                                           //needed for desired error rate. Test if N >= 3 is
return Witness(RandInt(2, N-2), N-1, N) == 1;    //prime using one value of A
}

int main()
{
srand(time(NULL));                //define random numbers
for(int i = 101; i<200; i+= 2)            //check primality of a few numbers
if(CheckPrime(i))                //if number is prime print it
cout <<  i << ” is a prime number” << endl;
}