This function applies the Fermat primality test on an integer, with a customizable amount of runs. The Fermat Test is based on Fermat's little theorem, which states that if $p$ is prime and $a$ is an integer, the following statement holds true:
$$a^p\equiv a\mod p.$$
If $a$ is not divisible by $p$, this is equivalent to
$$a^{p-1}\equiv 1\mod p.$$
By reversing this statement, it follows that if
$$a^{p-1}\not\equiv 1~\mod p$$
holds, $p$ can't be prime. By repeatedly checking this statement for different values of $a$, with $a\in[2,p-2]$, one can determine, whether $p$ is "likely prime" or not.
It should be noted that there exist numbers $n$, with $a^{n-1}\equiv 1\mod n$ for all $a<n$ with $\gcd(a,n)=1$. These numbers are called Carmichael numbers and they can only potentially be detected as non-prime by integers $a$ with $\gcd(a,n)\not=1$.
The implementation here uses the fundamental long long data type and a modular exponentiation function. Thus, it only works for numbers up to $\sqrt{2^{63}-1}\approx 3\cdot 10^9$. To avoid this restriction, a data type of unrestricted size could be used.
Download(zip).
//Knowledgedump.org - Function for modular exponentiation by squaring, to reduce integer overflows.
#ifndef MODULAR_EXP_H //Include guard
#define MODULAR_EXP_H
#include <iostream>
//Forward declaration
//Function for checking input arguments.
long long modular_exp(long long base, long long exp, long long mod);
//Recursive function for calculation.
long long mod_exp(long long base, long long exp, long long mod);
//Function only allows exponent bigger or equal 0, since integer types are used.
//Checks input, then calls the actually calculating function.
long long modular_exp(long long base, long long exp, long long mod) {
if (base == 0) { return 0; }
if (base == 1) { return 1; }
if (exp == 0) { return 1; }
if (exp < 0) {
std::cout << "Invalid exponent input." << std::endl;
return 0;
}
if ((mod == 0) || (mod == 1) || (mod == -1)) {
std::cout << "Invalid modulo input." << std::endl;
return 0;
}
return mod_exp(base, exp, mod);
}
//Since function is called recursively, the actual calculation process is done in a separate function,
//to avoid redundant base, exp and mod checks.
long long mod_exp(long long base, long long exp, long long mod) {
if (exp == 1) { return base % mod; }
if (exp > 0) {
if (exp % 2 == 0) {
long long mid = mod_exp(base, exp / 2, mod);
return (mid * mid) % mod;
}
else {
return (base * mod_exp(base, exp - 1, mod)) % mod;
}
}
}
#endif //Include guard
//Knowledgedump.org - Function that runs a Fermat Test, in order to detect, whether a number is prime.
//Works for numbers up to about 3*10^9.
#ifndef MEASURING_ERROR_H //Include guard
#define MEASURING_ERROR_H
#include <iostream>
#include <cstdlib> //For RNG functions rand and srand
#include <ctime> //For RNG seed time(0)
#include "modular_exp.h" //For pow
//n is the number to check for being prime and k is the amount of runs.
//long long chosen for allowance of larger numbers, and compatibility with modular_exp().
bool fermat_test(long long n, long long k) {
if (n <= 1) {
std::cout << "Number too small." << std::endl;
return false;
}
if ((n == 2) || (n == 3)) {
std::cout << n << " is prime." << std::endl;
return true;
}
if (k <= 0) {
std::cout << "Invalid amount of runs." << std::endl;
return false;
}
//Pick "random" seed.
std::srand(time(0)); //time(0) is the time difference between current time and 00:00, 01-01-1970 UTC in seconds.
long long a = 1;
while (k > 0) {
a = (std::rand() % (n - 3)) + 2; //a is random number between 2 and max(n-2, RAND_MAX+2).
if (modular_exp(a, n-1, n) != (1 % n)) { //modular_exp() instead of pow(), to avoid integer overflow.
std::cout << n << " is not prime." << std::endl;
return false;
}
--k;
}
std::cout << n << " is likely prime." << std::endl;
return true;
}
#endif //Include guard
//Knowledgedump.org - Example file.
#include "fermat_test.h"
int main() {
long long i = 1;
long long j = 1;
//Run fermat_test in a loop, terminate program with 0. For valid results, number has to be below 3*10^9.
while (i > 0) {
std::cout << "Enter number to apply fermat primality test on. 0 for program termination." << std::endl;
std::cin >> i;
if (i != 0) {
std::cout << "Enter number of runs." << std::endl;
std::cin >> j;
std::cout << "Result: ";
fermat_test(i,j);
}
}
return 0;
}
Enter number to apply fermat primality test on. 0 for program termination. 655360001 Enter number of runs. 1000 Result: 655360001 is likely prime. Enter number to apply fermat primality test on. 0 for program termination. 655360003 Enter number of runs. 1000 Result: 655360003 is not prime. Enter number to apply fermat primality test on. 0 for program termination. 0