# Deterministic Miller-Rabin implementation

I'm trying to implement a primality check function with a deterministic Miller-Rabin algorithm but the results are not always correct: when checking first 1,000,000 numbers it only founds 78,495 instead of 78,498.

This is obtained using [2, 7, 61] as a base which, according to wikipedia, should always be correct for values up to 4,759,123,141.
The interesting thing is that the 3 missing primes are exactly the ones componing the base (2, 7 and 61).

Why is this happening? The code I'm using is the following:

``````T modular_power(T base, T exponent, T modulo) {
base %= modulo;
T result = 1;

while (exponent > 0) {
if (exponent % 2 == 1)
result = (result * base) % modulo;
base = (base * base) % modulo;
exponent /= 2;
}

return result;
}

bool miller_rabin(const T& n, const vector<T>& witnesses) {
unsigned int s = 0;
T d = n - 1;
while (d % 2 == 0) {
s++;
d /= 2;
}

for (const auto& a : witnesses) {
if (modular_power<T>(a, d, n) == 1)
continue;

bool composite = true;
for (unsigned int r = 0; r < s; r++) {
if (modular_power<T>(a, (T) pow(2, r) * d, n) == n - 1) {
composite = false;
break;
}
}

if (composite)
return false;
}

return true;
}

bool is_prime(const T& n) {
if (n < 4759123141)
return miller_rabin(n, {2, 7, 61});
return false; // will use different base
}
``````