In-depth introduction to prime number Algorithms

// For the missing code, see </P> <p> # include <stdlib. h> </P> <p> static int CMP (const int * P, const int * q) <br/>{< br/> return (* P) -(* q); <br/>}</P> <p> bool isprime (int n) <br/>{< br/> If (n <2) return false; <br/> If (n = 2) return true; <br/> If (N % 2 = 0) return false; </P> <p> If (n> = 67 & n <= primes [nelems (primes)-1]) <br/>{< br/> return null! = <Br/> bsearch (& N, primes, nelems (primes), sizeof (N), CMP ); <br/>}< br/> else <br/> {<br/> for (INT I = 1; primes [I] * primes [I] <= N; + + I) <br/> If (N % primes [I] = 0) return false; <br/> return true; <br/>}< br/>

NOTE: If there is no special explanation, we will discuss the time complexity when n is a prime number.

1. Judge by concept:

If a positive integer has only two factors, 1 and P, P is called a prime number.

Bool isprime (int n) <br/>{< br/> If (n <2) return false; </P> <p> for (INT I = 2; I <n; ++ I) <br/> If (N % I = 0) return false; </P> <p> return true; <br/>}</P> <p>

Time complexity O (n ).

2. Improve and remove the even number judgment

Bool isprime (int n) <br/>{< br/> If (n <2) return false; <br/> If (n = 2) return true; </P> <p> for (INT I = 3; I <n; I + = 2) <br/> If (N % I = 0) return false; </P> <p> return true; <br/>}</P> <p>

Time complexity O (n/2), speed doubled.

3. Further reduce the scope of judgment

Theorem: If n is not a prime number, N has a factor D that satisfies 1 <D <= SQRT (n.
Proof: If n is not a prime number, the defined N has a factor D that satisfies 1 <D <n.
If D is greater than SQRT (N), N/D is a factor that satisfies 1 <n/d <= SQRT (n.

Bool isprime (int n) <br/>{< br/> If (n <2) return false; <br/> If (n = 2) return true; </P> <p> for (INT I = 3; I * I <= N; I + = 2) <br/> If (N % I = 0) return false; </P> <p> return true; <br/>}</P> <p>

You can also

Bool isprim (int A) <br/>{< br/> int divisor = 3; <br/> int Limit =; <br/> if (a % 2 = 0) <br/> return false; <br/> while (Limit> divisor) <br/>{< br/> if (a % divisor = 0) <br/> return false; <br/> Limit = A/divisor; <br/> divisor + = 2; <br/>}< br/> return true; <br/>}

Time complexity O (SQRT (N)/2), speed increase O (n-SQRT (N)/2 ).

4. Duplicate judgment in elimination factors.
Example: 11% 3! = 0 can be determined 11% (3 * I )! = 0.

Theorem: If n is not a prime number, N has a "prime number" factor D that satisfies 1 <D <= SQRT (n.
Proof: i1. if n is not a prime number, then n has a factor D that satisfies 1 <D <= SQRT (n.
I2. if D is a prime number, the theorem proves that the algorithm ends.
I3. make n = D and go to step i1.

Since it is impossible to infinitely break down n factors, the algorithm proved above will eventually stop.

// Primes [I] is an ascending sequence of prime numbers: 2, 3, 5, 7 ,... <br/> // more accurately, the PRIMES [I] sequence contains 1-> SQRT (N) all prime numbers in the range </P> <p> bool isprime (INT primes [], int N) <br/>{< br/> If (n <2) return false; </P> <p> for (INT I = 0; primes [I] * primes [I] <= N; ++ I) <br/> If (N % primes [I] = 0) return false; </P> <p> return true; <br/>}</P> <p>

Assume that the number of prime numbers in the N range is Pi (N), then the time complexity O (PI (SQRT (n ))).

The PI (x) function satisfies the prime number theorem: ln (x)-3/2 <X/PI (x) <Ln (x)-1/2, when x> = 67.

Therefore, O (PI (SQRT (N) can be expressed as O (SQRT (X)/(Ln (SQRT (x)-3/2 )),

O (SQRT (X)/(Ln (SQRT (x)-3/2) is also the spatial complexity of this algorithm.

5. Construct the PRIMES [I]: 2, 3, 5, 7 ,...

We know from the 4 algorithm that it is very efficient to judge whether N is a prime number when the prime number sequence has been constructed;

However, when constructing the prime number sequence itself, does it achieve the best efficiency?

In fact, this is acceptable! -- We can fully utilize the constructed prime number sequence during construction!

Let us assume that we already have a prime number sequence: P1, P2,... PN

To determine whether Pn + 1 is a prime number, all prime number sequences in the range of (1, SQRT (Pn + 1)] are required,

This prime number sequence is already included as a subset of P1, P2,... PN!

// Construct the PRIMES sequence [] </P> <p> void makeprimes (INT primes [], int num) <br/>{< br/> int I, j, CNT; </P> <p> primes [0] = 2; <br/> primes [1] = 3; </P> <p> for (I = 5, CNT = 2; CNT <num; I + = 2) <br/>{< br/> int flag = true; <br/> for (j = 1; primes [J] * primes [J] <= I; ++ J) <br/> {<br/> if (I % primes [J] = 0) <br/>{< br/> flag = false; break; <br/>}< br/> If (FLAG) primes [CNT ++] = I; <br/>}</P> <p>

The time complexity of makeprimes is complex, and it is called only once during initialization. in a certain range of applications, we can consider that makeprimes requires constant time. in the subsequent discussions, we will discuss a better makeprimes Method for computers.

Appendix: Deletion of prime numbers

[Theorem] If the multiples of all prime numbers smaller than prime p have been removed from isprim and p * P> N, all the remaining numbers are prime numbers.

 Void creatprim () <br/>{< br/> for (INT I = 2; I <= N; I ++) <br/> isprim [I] = I % 2; <br/> isprim [2] = 1; <br/> for (INT I = 3; I * I <n; I + = 2) <br/>{< br/> If (isprim [I]) <br/> for (Int J = I; I * j <= N; j + = 2) <br/> isprim [I * j] = 0; <br/>}</P> <p>
6. Better use of computer resources...

In the current mainstream PC, the size of an integer is 2 ^ 32. if you need to determine whether the number of 2 ^ 32 is a prime number, you may need to test [2, all prime numbers in the range of 2 ^ 16] (2 ^ 16 = SQRT (2 ^ 32 )). from the Prime Number Theorem mentioned in 4, we can roughly determine the number of prime numbers in the range [2, 2 ^ 16. because 2 ^ 16/(Ln (2 ^ 16)-1/2) = 6138, 2 ^ 16/(Ln (2 ^ 16)-3/2) = 6834, we can estimate the number of prime numbers in the range of [2, 2 ^ 16] By 6138 <Pi (2 ^ 16) <6834.

After counting the prime numbers in the range of [2, 2 ^ 16], we found that there are only 6542 prime numbers:

P_6542: 65521,655 21 ^ 2 = 4293001441 <2 ^ 32, (2 ^ 32 = 4294967296)
P_6543: 65537,655 37 ^ 2 = 4295098369> 2 ^ 32, (2 ^ 32 = 4294967296)

In the actual operation, unsigned long x = 4295098369; will overflow, which is 131073. in the program, I calculated the result using the double type.

From the analysis, we can see that we only need to buffer 6543 prime numbers, and we can use the algorithm in 4 to efficiently judge the prime numbers in such a large range as [2, 2 ^ 32! (The size of the original 2 ^ 32 problem has now been reduced to 6543 !)

Although there is no problem with the current computer processing the 6542 prime numbers in the range of [2, 2 ^ 16], although makeprimes only needs to be run once, but we still need to consider whether it is possible to be improved ?!

I want to use makeprimes as a static initialization implementation. In C ++, I can simulate a similar implementation of static initialization in Java:

# Define nelems (x) (sizeof (x)/(sizeof (x) [0])

Static int primes [6542 + 1];
Static struct _ init {_ Init () {makeprimes (primes, nelems (primes) ;}}_ Init;

In this way, you can automatically initialize the prime number sequence with makeprimes when the program starts. However, my current idea is: why cannot we call the makeprimes function during compilation? Yes !!! The Code is as follows:

// This code can be directly generated by the Program </P> <p> const static int primes [] = <br/> {<br/>, 97,101,103, 107,109,113,127,131,137,139,149,151,157,163,167,173,179,181,191,193,197,199,211, 223,227,229,233,239,241,251,257,263,269,271,277,281,283,293,307,311,313,317,331, 337,347,349,353,359,367,373,379,383,389,397,401,409,419,421,431,433,439,443,449, <br/>, <br/> 457,461,463,467,479,487,491,499,503,509,521,523,541,547,557,563,569,571,577,587, <br/> 593,599,601,607,613,617,619,631,641,643,647,653,659,661,673,677,683,691,701,709, <br/> 719,727,733,739,743,751,757,761,769,773,787,797,809,811,821,823,827,829,839,853, <br/> 857,859,863,877,881,883,887,907,911,919,929,937,941,947,953,967,971,977,983,991, <br/>... <br/> 65521,655 37 <br/>}; </P> <p>

A bit incredible: the time complexity originally required by makeprimes has now actually changed to O (1! (I think o (0) may be more suitable !)

7. Binary Search

Now, We cache the list of first SQRT (2 ^ 32)/(Ln (SQRT (2 ^ 32)-3/2) prime numbers.

When a prime number is used, only Pi (SQRT (2 ^ 32) queries are required at most (the exact value is 6543 times). Is there any other way to determine the number?

When the prime number is small (not greater than 2 ^ 16), can I directly query the cached Prime Number list?

The answer is yes! Because primes is an ordered sequence, when the prime number is smaller than 2 ^ 16, we can directly

The binary method is used to obtain from PRIMES (if the query fails, it is not a prime number ).

 

Time Complexity:

If (n <= primes [nelems (primes)-1] & n> = 67): O (log2 (nelems (primes) <13;
If (n> primes [nelems (primes)-1]): O (PI (SQRT (N) <= nelems (primes ).

8. Prime Number Theorem + 2 Division Lookup

In 9, we use a two-way lookup method to judge the number of small values equal to primes [nelems (primes)-1. we have previously determined 13 times for the 6453 prime numbers buffered by 2 ^ 32 (log2 (1024*8) = 13 ). for a small prime number (in fact, it is only the number in the range of 2 ^ 16), 13 comparisons are completely acceptable. however, according to the prime number theorem: ln (x)-3/2 <X/PI (x) <Ln (x)-1/2, when x> = 67, we can still narrow down the search range smaller than 2 ^ 32 (from 0 to nelems (primes)-1 ). we need to solve the problem n <= primes [nelems (primes)-1 ):

If n is a prime number, where is the possible range of N in the prime number sequence?

---- (N/(Ln (N)-1/2), N/(Ln (N)-3/2), that is, the prime number theorem!

The above code is modified as follows:

Bool isprime (int n) <br/>{< br/> If (n <2) return false; <br/> If (n = 2) return true; <br/> If (N % 2 = 0) return false; </P> <p> int HI = (INT) Ceil (N/(Ln (N) -3/2); </P> <p> If (n> = 67 & Hi <nelems (primes )) <br/>{< br/> int Lo = (INT) floor (N/(Ln (N)-1/2); </P> <p> return null! = <Br/> bsearch (& N, PRIMES + lo, Hi-lo, sizeof (N), CMP ); <br/>}< br/> else <br/> {<br/> for (INT I = 1; primes [I] * primes [I] <= N; + + I) <br/> If (N % primes [I] = 0) return false; <br/> return true; <br/>}</P> <p>

Time Complexity:

If (n <= primes [nelems (primes)-1] & n> = 67): O (log2 (Hi-Lo) <???;
If (n> primes [nelems (primes)-1]): O (PI (SQRT (N) <= nelems (primes ).

9. Package the package into a prime number Library (giving all the code)

So far, I have provided all the methods I know for improvement (if someone has a better algorithm, thank you for telling me ). here, we need to emphasize that the prime number method discussed here is for the number in the range of 0-2 ^ 32. As for the number that looks for hundreds of bits, we will not discuss the scope here, it should be purely mathematical content.

The code is saved in two files: Prime. H, prime. cpp.

// File: Prime. h </P> <p> # ifndef prime_h_2006_10_27 _ <br/> # define prime_h_2006_10_27 _ </P> <p> extern int prime_max (void ); // size of the prime number sequence <br/> extern int prime_get (int I); // returns the I-th prime number, 0 <= I <prime_max </P> <p> extern bool prime_test (INT N); // test whether it is a prime number, 1 <= n <int_max </P> <p> # endif </P> <p> ////////////////// /// // </P> <p> // file: prime. CPP </P> <p> # include <assert. h> <br/> # include <Limits. h> <br/> # include <math. h> <br/> # include <stdlib. h> </P> <p> # include "Prime. H "</P> <p> // calculate the number of elements in the array </P> <p> # define nelems (x) (sizeof (x )) /(sizeof (x) [0]) </P> <p> // a sequence of prime numbers. It stores at least the first 6543 prime numbers! <Br/> // before I personally get used to buffering (1024*8) -) </P> <p> static const int primes [] = <br/> {<br/>, 97,101,103, 107,109,113,127,131,137,139,149,151,157,163,167,173,179,181,191,193,197,199,211, 223,227,229,233,239,241,251,257,263,269,271,277,281,283,293,307,311,313,317,331, 337,347,349,353,359,367,373,379,383,389,397,401,409,419,421,431,433, <br/>, 439,443,449, <br/> 457,461,463,467,479,487,491,499,503,509,521,523,541,547,557,563,569,571,577,587, <br/> 593,599,601,607,613,617,619,631,641,643,647,653,659,661,673,677,683,691,701,709, <br/> 719,727,733,739,743,751,757,761,769,773,787,797,809,811,821,823,827,829,839,853, <br/> 857,859,863,877,881,883,887,907,911,919,929,937,941,947,953,967,971,977,983,991, <br/>... <br/> 65521,655 37 <br/>}; </P> <p> // Bsearch comparison function </P> <p> static int CMP (const void * P, const void * q) <br/>{< br/> return (* (int *) P)-(* (int *) q ); <br/>}</P> <p> // Number of buffered prime numbers </P> <p> int prime_max () <br/>{< br/> return nelems (primes ); <br/>}</P> <p> // returns the I-th prime number </P> <p> int prime_get (int I) <br/> {<br/> assert (I> = 0 & I <nelems (primes); <br/> return primes [I]; <br/>}</P> <p> // test whether N is a prime number </P> <p> bool prime_test (int n) <br/>{< br/> assert (n> 0); </P> <p> // separate determination of even numbers </P> <p> If (n <2) return false; <br/> If (n = 2) return true; <br/> If (! (N & 1) return false; </P> <p> // if n is a prime number, it is prior to the sequence Hi position </P> <p> int Lo, hi = (INT) Ceil (N/(log (N)-3/2. 0); </P> <p> If (Hi <nelems (primes )) <br/> {<br/> // determine the range of the binary method search <br/> // The Prime Number Theorem is satisfied only when n> = 67 is used. </P> <p> if (n> = 67) lo = (INT) floor (N/(log (N)-1/2. 0); <br/> else {Lo = 0; Hi = 19 ;} </P> <p> // if the search is successful, the prime number is displayed. </P> <p> return null! = <Br/> bsearch (& N, PRIMES + lo, Hi-lo, sizeof (N), CMP ); <br/>}< br/> else <br/>{< br/> // events that are not within the range of the stored Prime Number Sequence </P> <p> for (int I = 1; primes [I] * primes [I] <= N; ++ I) <br/> If (N % primes [I] = 0) return false; </P> <p> return true; <br/>}</P> <p>
10. Review and promotion

The discussion on prime numbers has come to an end. looking back at our previous solutions, we will find that it is difficult to design good algorithms without the basic knowledge of mathematics. However, if we only consider the mathematical principles, the same problem also arises when computer essential features are neglected. A common example is to calculate the Fibonacci series. of course there are many methods, but there is no need to implement them in the current computer!

Because the Fibonacci series itself increase exponentially, 32-bit signed integers can represent only the first 46:

Static const int maid [] = <br/> {<br/> 144,233,377,610,987,159, 17711,28657, 462.16,75025, 121393,196418, <br/> large, 832040,1346269, large, large, <br/> large, 39088169,63245986, large, 267914296, <br/> large, small, -1323752223 </P> <p> // Note: The f47 series has exceeded !!! <Br/>}; </P> <p>

Therefore, I only need to save the first 46 pieces of Fiber-ACCI to the array! For example, F (int I) {return Fibonacci (I);} is very simple and efficient. For example, find a factorial N !, Although there are also a lot of mathematical descriptions, the computer generally cannot express it. We can simply put usable computing into the memory.

In short, many things are good, but do not be bound by them!