Today, I will talk about a more advanced algorithm with a simple purpose: finding prime numbers. However, we will be using a new algorithm called the Miller-Rabin algorithm, which is a probabilistic algorithm based on the Fermat's little theorem. It may sound mysterious, but I have just recently come across this algorithm myself. It is a purely mathematical solution, so if you don't understand it, it might be a good idea to learn some mathematics. Let's continue.
First, let's understand the basic mathematical knowledge, Fermat's little theorem:
If n is a prime number, then for all integers a such that 1≤a≤n-1, we have a^(n-1)mod n=1.
The author Fermat was a great mathematician, and he lived in a time when computers did not exist. However, today we will use this theorem to design an algorithm.
By analyzing this theorem, we can conclude that if a number is a prime number, then it must satisfy the condition that for any integer a in the range (1, n-1), a^(n-1)mod n=1. Don't understand? Let's try the contrapositive: if there exists an integer a in the range (1, n-1) such that a^(n-1)mod n=1 does not hold, then the number is not a prime number. I believe you understand now. To determine whether a number is prime, we just need to randomly generate a series of numbers a. If these numbers a satisfy Fermat's little theorem for the number N, then we can consider it as a prime number. Of course, there is a prerequisite: this deduction can only be used in the field of computers. It's similar to the optimization of an algorithm problem using the Goldbach conjecture that I mentioned last time. It can only be used within the range of numbers that computers can calculate. Please do not use it in other scientific fields.
The algorithm is as follows. Originally, it was a C++ algorithm I found on the internet, but I modified it to C.
/*
- Application of Fermat's little theorem for finding prime numbers
- Miller-Rabin algorithm
- 2008 12 27
*/
#include"stdio.h"
#include"stdlib.h"
#include"math.h"
#include"time.h"
int Btest(int a,int n)
{
int s = 0;
int t = n-1;
int i , j , k;
int x = 1;
int y;
i = 1;
do{
s++;
t = t / 2;
}while((t%2)!=1);
while(i <= t){
x = ( x * a ) % n;
i++;
}
if((x == 1) || ( x == n-1))
return 1;
for(j = 1 ; j <= s - 1 ; j++){
y = 1;
for(k = 1 ; k <= j ; k++)
{
y = 2 * y;
}
i = 1;
x = 1;
while(i <= ( y * t)){
x = ( x * a) % n;
i++;
}
if(x == n - 1)
return 1;
}
return 0;
}
int MillRab(int n)
{
int a;
/* Use time as a random seed to generate random numbers */
a=random((unsigned)time(0))%(n-3)+2;
return Btest(a,n);
}
int MillerRabin(int n,int k)
{
int i;
for(i=1;i<=k;i++)
{
if(MillRab(n)==0)return 0;
}
return 1;
}
int main(){/* Application of Fermat's little theorem /
int i;
int n=10000;/ Define the testing range /
int count=0;
printf("2 3 ");/ Print 2 and 3 first /
for(i=5;i<=n;){/ Loop starting from 5 /
if(MillerRabin(i , (int)log10(i))){/ Satisfies the condition? */
printf("%d ",i);
count++;
}
i=i+2;
}
printf("nthere %d prime(s)n",count);
return 0;
}