A small implementation of the Miller-Rabin primality test. This version is deterministic and correct for numbers up to 2,152,302,898,747.

MillerRabin.cpp

typedef long long tint;

tint powmod(tint a, tint k, tint n) {
  tint b = 1;
  while (k) {
    if (k & 1) {
      --k;
      b = (b * a) % n;
    } else {
      k >>= 1;
      a = (a * a) % n;
    }
  }
  return b;
}

bool miller_rabin_round(tint n, tint d, int s, tint a) {
  tint x = powmod(a, d, n);
  if (a < 2 || a > n - 2) return true;
  if (x == 1 || x == n - 1) return true;
  for (int i = 1; i < s; ++i) {
    x = (x * x) % n;
    if (x == 1) return false;
    if (x == n - 1) return true;
  }
  return false;
}

bool miller_rabin(tint n) {
  if (n <= 2) return n == 2;
  tint nm = n - 1;
  int s = __builtin_popcount(nm  ^ (nm - 1)) - 1;
  tint d = nm >> s;
  return miller_rabin_round(n, d, s, 2) &&
         miller_rabin_round(n, d, s, 7) &&
         miller_rabin_round(n, d, s, 61);
}