徹底理解線性篩選法


問題:求所有小於等於n(n比較大)的所有素數

首先可能最容易想到的是寫一個函數來判斷它是不是素數,但是對於求比n小的整數就顯得時間復雜度太高了,一般解這種問題會采用篩選法...

埃氏篩選法

思想是,使用一個位數組is_prime保存每一個數是否是素數,然后每次找到一個素數x,就把這個素數的i倍,i滿足(x * i <= n)踢掉,即將is_prime設置為true(默認是false)。這樣找到一個數后,再找的下一個is_prime = false的數,一定是素數(因為如果不是,它就會被某個比它小的數因子踢出掉)

線性篩選法

埃氏篩選法的效率已經不低了,但是會有很多的重復計算,即一個數會被多次踢出,浪費了計算時間。對於任何一個合數z,一定存在質數x,使得:

z = x * y (y 不為 1)

但是這樣的x不唯一,比如6=23=32,這里x = 2或者3。所以我們現在要找一個辦法,讓每一個合數只能被踢出一次而不是多次。這時注意到x的最小值是一定的,所以我們可以用每一個合數的最小質因子來踢出。

即:

10由2踢出
9由3踢出
77由7踢出

這樣就實現了一個合數被一個最小質數踢出,大大加快了運算效率。

那么在程序實現時怎么做呢?
我們可以在遍歷2-n時,用當前已經找到的每一個質數去乘以這個當前數i,在這些質數眼里(不是當前遍歷的數!!),這樣就相當於:

質數2分別與2,3,4,5,6,7,8,9...,踢出了(4,6,8,10,12,14,16,18)
質數3分別與3,4,5,6,7,8,9...踢出了(9,12....)
質數5分別與5,6,7,8,9,10,11,12...,踢出了(25,30...)
質數7分別與3,4,5,6,7,8,9...踢出了(21,28....)

這樣看起來還是有不少重復啊,不要急,目前這種方式相當於埃氏篩選法,這里最核心的一步是,如果當前數i的因子已經有了某個質數的話,那么當前數i就不再繼續與下一個質數相乘了。

什么意思?
就是說,當i = 6時(遍歷到6這個合數時),我先把2 * 6 = 12踢掉,當我想繼續再想把3(第二個質數) * 6 = 18頁踢掉時,這時候等一等!由於:

6 % 2 == 0

這個時候就不要繼續往下走了,因為將要計算的下一個要踢出的數的質最小因子一定是剛剛那個質數(2)!,為什么?因為當前數都可以整除剛剛那個質數了!

這樣就保證了一個合數只能被一次踢出,比如:

100 一定是在i = 50被踢掉的
70 一定是在i = 35被踢掉的
45 一定是在i = 15被踢掉的

證明:
首先證明合數t一定會被踢掉:

對於t,一定存在最小質數a,所以當遍歷到i = t / a 時,踢掉t/a * a = t

然后證明t只會被踢一次

使用反證法,設i = x,i = y時踢出t,這里x < y,a、b是質數且
a是最小質因子(已證明這個一定存在),則

t = x * a
t = y * b

由於任意一個合數都可以表示為am*bn.....其中a,b,c...是質數,如45=3^2*5
所以當b!=a時,b一定大於a,這個時候循環已經超過a了(a是最小質因子),而45已經含有a這個因子,所以有t % a == 0得,此時循環已經break,因此不會重復計算任意合數。

代碼

埃氏篩選法:

#include <iostream>
#include <vector>
using namespace std;

vector<int> primes;

bool *is_prime;

void get_all_prime(int n){
	is_prime = new bool[n + 1];
	
	for(int i = 2; i <= n; i++){
		if(!is_prime[i]){
			primes.push_back(i);
			
			for(int j = i * 2; j <= n; j = j + i){
				is_prime[j] = true;
			}
		}
	}
}

int main(){
	get_all_prime(1000000000);
	
	//int len = primes.size();
	//for(int i = 0; i < len; i++){
	//	cout << primes[i] << endl;
	//}
	
	cout << primes.size();
	
	return 0;
} 

計算[2,1000000000]素數花了26s
不使用vector花了24s

線性篩選法(使用vector)

#include <iostream>
#include <vector>
using namespace std;

vector<int> primes;

bool *is_prime;

void get_all_prime(int n){
	is_prime = new bool[n + 1];
	
	for(int i = 2; i <= n; i++){
		if(!is_prime[i]){
			primes.push_back(i);
		}
		
		int size = primes.size();
		for(int j = 0; j < size && primes[j] * i <= n; j++){
			is_prime[primes[j] * i] = true;
			if(i % primes[j] == 0)
				break;
		}
	}
}

int main(){
	get_all_prime(1000000000);
	
	//int len = primes.size();
	//for(int i = 0; i < len; i++){
	//	cout << primes[i] << endl;
	//}
	
	cout << primes.size();
	
	return 0;
} 

計算[2,1000000000]素數花了19s

線性篩選法(不使用vector)

#include <iostream>
#include <vector>

#define MAX_CAP 100000000
using namespace std;

bool *is_prime;
int data[MAX_CAP];
int size = 0;

void get_all_prime(int n){
	if(n == 2){
		data[0] = 2;
		size = 1;
		return;
	}
	
	if(n == 3){
		data[0] = 2;
		data[1] = 3;
		size = 2;
		return;
	}
	
	is_prime = new bool[n + 1];
	
	for(int i = 2; i <= n; i++){
		if(!is_prime[i]){
			data[size++] = i;
		}
		
		for(int j = 0; j < size && data[j] * i <= n; j++){
			is_prime[data[j] * i] = true;
			if(i % data[j] == 0)
				break;
		}
	}
}

int main(){
	get_all_prime(1000000000);
	
	cout << size;
	
	return 0;
} 

計算[2,1000000000]素數花了9s。。。


免責聲明!

本站轉載的文章為個人學習借鑒使用,本站對版權不負任何法律責任。如果侵犯了您的隱私權益,請聯系本站郵箱yoyou2525@163.com刪除。



 
粵ICP備18138465號   © 2018-2025 CODEPRJ.COM