# | Likes | Tech tags | Title | Creator | Created date |
---|---|---|---|---|---|
1 | 0 |
2022-11-26 19:19
|
A simple implementation of the Miller-Rabin primality test algorithm for 64 bit integers. Runs in \(O(k \log^2 n)\) thanks to the fact that modular addition is constant for 64-bit integers (\(k\) is the number of tests, \(n\) is the integer value).
Caveat: The value of \(2n\) must not overflow a 64-bit integer. For larger integers, you need to use some big-integer class with addition and modulo operations implemented
#include<random>
long long modmult(long long a, long long b, long long mod) { // Compute a*b % mod
long long result = 0;
while(b > 0) {
if(b%2 == 1) {
result = (result+a) % mod;
}
a = (a+a) % mod;
b /= 2;
}
return result;
}
long long modpow(long long a, long long b, long long mod) { // Compute a^b % mod
long long result = 1;
while(b > 0) {
if(b%2 == 1) {
result = modmult(result, a, mod);
}
a = modmult(a, a, mod);
b /= 2;
}
return result;
}
// Miller–Rabin primality test
bool prime_test(long long n, int num_tests = 20) { // Assumes n is an odd integer larger than 3
static std::mt19937_64 randgen(0); // Change seed if needed
long long d = n-1;
long long s = 0;
while(d%2 == 0) {
s++;
d /= 2;
}
for(int test = 0; test < num_tests; test++) {
long long a = randgen()%(n-3) + 2; // Between 2 and v-2
long long x = modpow(a, d, n);
for(int i=0; i<s; i++) {
long long y = modmult(x, x, n);
if(y == 1 && x != 1 && x != n-1) { // Nontrivial square root of 1 modulo n
return false; // (x+1)(x-1) divisible by n, meaning gcd(x+1, n) is a factor of n, negating primality
}
x = y;
}
if(x != 1) {
return false;
}
}
return true; // Number is prime with likelihood of (1/4)^num_tests
}
int main() {
prime_test(5127821565631733);
}
objects | ||
std::mt19937_64 |
en.cppreference.com | cplusplus.com |
functions | ||
std::mersenne_twister_engine::mersenne_twister_engine |
en.cppreference.com | cplusplus.com |
std::mersenne_twister_engine::operator() |
en.cppreference.com | cplusplus.com |
Implement the Miller-Rabin primality test algorithm. Use your implementation to test the primality of some number.