原理：

　　素數，指在一個大於1的自然數中，除了1和此整數自身外，不能被其他自然數整除的數。在加密應用中起重要的位置，比如廣為人知的RSA演算法中，就是基於大整數的因式分解難題，尋找兩個超大的素數然後相乘作為金鑰的。一個比較常見的求素數的辦法是埃拉託斯特尼篩法(the Sieve of Eratosthenes) ，說簡單一點就是畫表格，然後刪表格，如圖所示：

　　從2開始依次往後面數，如果當前數字一個素數，那麼就將所有其倍數的數從表中刪除或者標記，然後最終得到所有的素數。

有一個優化：

標記2和3的倍數的時候，6被標記了兩次。所以從i的平方開始標記，減少很多時間。

比如3的倍數從9開始標記，而不是6，並且每次加6。

除了2以外，所有素數都是奇數。奇數的平方還是奇數，如果再加上奇數就變成了偶數一定不會是素數，所以加偶數(2倍素數)。

預先處理了所有偶數。

注意：1既不是素數也不是合數，這裡沒有處理1。

#! prime.py
import time

def primes(n):
  P = []
  f = []
  for i in range(n+1):
    if i > 2 and i%2 == 0:
      f.append(1)
    else:
      f.append(0)

  i = 3
  while i*i <= n:
    if f[i] == 0:
      j = i*i
      while j <= n:
        f[j] = 1
        j += i+i
    i += 2

  P.append(2)
  for i in range(3,n,2):
    if f[i] == 0:
      P.append(i)

  return P

def isPrime(n):
  if n > 2 and n%2 == 0:
    return 0

  i = 3
  while i*i <= n:
    if n%i == 0:
      return 0
    i += 2

  return 1

def primeCnt(n):
  cnt = 0
  for i in range(2,n):
    if isPrime(i):
      cnt += 1
  return cnt

if __name__ == '__main__':
  start = time.clock()
  n = 10000000
  P = primes(n);
  print("There are %d primes less than %d"%(len(P),n))
  #for i in range(10):
  #  print(P[i])
  print("Time: %f"%(time.clock()-start))
  #for n in range(2,100000):
  #  if isPrime(n):
  #    print("%d is prime"%n)
    #print("%d is "%n + ("prime" if isPrime(n) else "not prime"))

  start = time.clock()
  n = 1000000
  print("There are %d primes less than %d"%(primeCnt(n),n))
  print("Time: %f"%(time.clock()-start)

 

用素數篩選法求1千萬以內的素數用了5.767s，

普通素數判斷法求1百萬以內的素數用了9.642s，

用C++素數篩選法求1億以內的素數用了0.948s，

用C++普通素數判斷法求1千萬以內的素數用了3.965s，

可見解釋語言確實比編譯語言慢很多。

附C++程式，用了位壓縮優化空間

#include <iostream>
#include <cstdio>
#include <algorithm>
using namespace std;
#define N 100000001

unsigned f[(N>>5)+5];
int p[5761456],m;
void init()
{
	int i,j;
	for(i=4;i<N;i+=2)
		f[i>>5]|=1<<(i&0x1F);
	p[m++]=2;
	for(i=3;i*i<N;i+=2)
		if(!(f[i>>5]&(1<<(i&0x1F))))
		{
			p[m++]=i;
			for(j=i*i;j<N;j+=i+i)
				f[j>>5]|=1<<(j&0x1F);
		}
	for(;i<N;i+=2)
		if(!(f[i>>5]&(1<<(i&0x1F))))
			p[m++]=i;
}
int is_prime(int n)
{
	int i;
	for(i=0;p[i]*p[i]<=n;i++)
		if(n%p[i]==0)
			return 0;
	return 1;
}
int isPrime(int n)
{
	if(n>2 && n%2==0)
		return 0;
	int i=3;
	while(i*i<=n)
	{
		if(n%i==0)
			return 0;
		i+=2;
	}
	return 1;
}
int main()
{
	int n=0,i;
	clock_t st=clock();
	init();
	/*for(i=2;i<10000000;i++)
		if(isPrime(i))
			n++;*/
	printf("%d %dms\n",m,clock()-st);
	/*while(~scanf("%d",&n),n)
	{
		i=lower_bound(p,p+m,n+1)-p;
		printf("%d\n",i);
	}*/
    return 0;
}