1. 程式人生 > >O(NloglogN)素數篩法與O(N)素數篩法的對比測試

O(NloglogN)素數篩法與O(N)素數篩法的對比測試

#include <ctime>
#include <cmath>
#include <iostream>
using namespace std;
const int mx = 1000000 + 1; ///在(1,mx)的範圍內尋找素數
const int sqrt_mx = (int)sqrt((double)mx);

bool vis[mx];
int prime[mx / 10]; ///在mx>65000時建議寫成 int prime[mx/10];

/*複雜度O(NloglogN)*/
void getPrime(int &cnt)
{
	int i, j;
	for (i = 2; i <= sqrt_mx; ++i)
		if (!vis[i])
		{
			prime[cnt++] = i;
			for (j = i * i; j < mx; j += i) vis[j] = true;
		}
	for (; i < mx; ++i) if (!vis[i]) prime[cnt++] = i;
}

/*複雜度O(N),但常數比較大
O(N)是因為每個合數至多被賦true值一次(每個合數n僅在i=n的最大非平凡因子(不為1和n的因子)時被賦值)
以12為例,在原來的演算法中,vis[12]在i=2和i=3時均被賦值為true
而在此演算法中,vis[12]只在i=6時被賦值為true(6為12的最大非平凡因子)
*/
void linear_getPrime(int &cnt)
{
	int i, j;
	for (i = 2; i < mx; ++i)
	{
		if (!vis[i]) prime[cnt++] = i;
		for (j = 0; j < cnt && i * prime[j] < mx; ++j)
		{
			vis[i * prime[j]] = true;
			if (i % prime[j] == 0) break;
		}
	}
}

/*
(測試資料因機器而異)

C++編譯器:
mx      getPrime()		linear_getPrime()		加速比
1e6		12.231ms		9.943ms					1.230
1e7     203.71ms		105.23ms				1.936
1e8		2319.3ms		1086.5ms				2.135
1e9		26762ms			11931ms					2.243

但是,在G++編譯器下:
mx      getPrime()		linear_getPrime()		加速比
1e6		4.37ms		    6.83ms					*0.640
3e6     20.51ms         20.19ms                 1.016
1e7     154ms		    76.8ms				    2.005
1e8		1885ms		    815.875ms				2.301
1e9		22258ms			9335ms					2.384

*/
void run()
{
	int cnt = 0;
	getPrime(cnt);
	//linear_getPrime(cnt);
}

int main()
{
	int loop = 100;
	double starttime = clock();
	//
	for (int i = 0; i < loop; ++i)
		run();
	//
	double endtime = clock();
	double usetime = (double)difftime(endtime, starttime) / loop;
	cout << "runtime:" << usetime << "ms" << endl;
	return 0;
}


結論:對於G++編譯器的使用者來說(福利),O(N)的演算法僅在mx較大時才有很明顯的優化效果

PS:還有一種線性篩法的實現,但測試中執行時間比上面那個多20%

const int mx = 100000000 + 1; ///在(1,mx)的範圍內尋找素數

int prime[mx / 10];
int nowPrimeNumber[mx];

void linear_getPrime(int &cnt)
{
	int i, j;
	for (i = 2; i < mx; ++i)
	{
		if (!nowPrimeNumber[i]) prime[++cnt] = i, nowPrimeNumber[i] = cnt;
		for (j = 1; j <= nowPrimeNumber[i] && i * prime[j] < mx; ++j)
			nowPrimeNumber[i * prime[j]] = j;
	}
}