线性筛素数详细整理

Posted 2020-12-05 lyfoi

如果你在1个月前让我判断素数，我一定会猛敲出以下代码：

bool check( int num )  
{  
     int tmp =sqrt( num);  
     for(int i= 2;i <=tmp; i++)  
        if(num %i== 0)  
          return 0 ;  
     return 1 ;  //实在是太慢了！
}  

$ $

下面给大家带来3种筛选素数和一种直接判断素数

$ $

---

$ $

什么是线性筛？

• 对于求多个质数时与其一个个判断不如用排除法，用空间换取大量时间。

$ $

---

$ $
$ $

一般筛法（埃拉托斯特尼筛法）：

基本思想

素数的倍数一定不是素数

实现方法

用一个长度为(N+1)的数组保存信息（(0)表示素数，(1)表示非素数），先假设所有的数都是素数（初始化为(0)），从第一个素数(2)开始，把(2)的倍数都标记为非素数（置为(1)），一直到大于(N)；然后进行下一趟，找到(2)后面的下一个素数(3)，进行同样的处理，直到最后，数组中依然为0的数即为素数

说明：整数(1)特殊处理即可。

举个例子

我们筛前(20)个数

首先初始为（(0)代表不是素数，(1)代表是素数）

(0) (1) (1) (1) (1) (1) (1) (1) (1) (1) (1) (1) (1) (1) (1) (1) (1) (1) (1) (1)

然后从(2)开始我们发现(2)被标记为素数，我们把(2)的倍数全部筛掉

变为：

(0) (1) (1) (0) (1) (0) (1) (0) (1) (0) (1) (0) (1) (0) (1) (0) (1) (0) (1) (0)

接着到(3)我们发现(3)仍然被标记，把(3)的倍数全部筛掉

变为：

(0) (1) (1) (0) (1) (0) (1) (0) (0) (0) (1) (0) (1) (0) (0) (0) (1) (0) (1) (0)

接着一直重复下去就得到了最后的素数表：

(0) (1) (1) (0) (1) (0) (1) (0) (0) (0) (1) (0) (1) (0) (0) (0) (1) (0) (1) (0)

(2) (3) (5) (7) (11) (13) (17) (19)

Code

const int MAXN = 1000000；  
void get_list()  
{  
    int i, j;  
    for (i=0; i<MAXN; i++) prime[i] = 1;  
    prime[0] = prime[1] = 0;  
    for (i=2; i<MAXN; i++)  
    {  
        if (!prime[i]) continue;  
        for (j=i*2; j<MAXN; j+=i) prime[ j ] = 0;  
    }  
}//调和级数证明可得复杂度为（nlogn）,所以不能称之为线性筛，但是它的实际运行速度也不是特别慢~~

$ $

---

$ $

$ $

下面我们来介绍一波真正的线性筛（欧拉筛法）：

我们发现在上面的筛法中有的数字是多个素数的倍数，也就是说它可能会被重复计算多次，比如说(6)同时是(2)与(3)的倍数，它在计算时就被访问了两次,这样会导致效率低下，所以在下面的算法中我们考虑如何优化这种情况。

基本思想

每一个合数可以被唯一地表示成它的一个最小质因子和另外一个数的乘积。

一个合数((x))与一个质数((y))的乘积可表示成一个更大的合数((Z))与一个更小的质数((a))的乘积，那样我们到每一个数，都处理一次，这样处理的次数是很少的，因此可以在线性时间内得到解。

实现方法

我们可以从2开始通过乘积筛掉所有的合数。将所有合数标记，保证不被重复筛除。

举个例子

仍然按上面的例子模拟（这里(0)为是素数，(1)为非素数，(p)为记录的素数表）：

初始：

(1) (0) (0) (0) (0) (0) (0) (0) (0) (0) (0) (0) (0) (0) (0) (0) (0) (0) (0) (0)

(p(empty))

然后到(2)的位置，把(2)放入素数表,做当前范围内可以筛掉的处理（具体是怎样的看代码叭）：

(1) (0) (0) (1) (0) (0) (0) (0) (0) (0) (0) (0) (0) (0) (0) (0) (0) (0) (0) (0)
(p) (2) 到(3)，把(3)放入素数表，继续处理

(1) (0) (0) (1) (0) (1) (0) (0) (1) (0) (0) (0) (0) (0) (0) (0) (0) (0) (0) (0)

(p) (2) (3) 然后到了(4)，它不是个素数，也处理一下

(1) (0) (0) (1) (0) (1) (0) (1) (1) (0) (0) (0) (0) (0) (0) (0) (0) (0) (0) (0)

(p) (2) (3) .......

然后一直搞下去，最后也能得到完整的素数表，这样虽然看起来复杂一些，但是实际上我们发现对于每个数的处理几乎是(O(1))的。

Code:

 void get_list(){
       for(int i=2;i<=maxn;i++){
             if(!is_not_pr[i]) prime[++tot]=i;
             for(int j=1;j<=tot&&i*prime[j]<=maxn;j++){
                   is_not_pr[i*prime[j]]=1;//合数标为1，同时，prime[j]是合数i*prime[j]的最小素因子
                   if(i%prime[j]==0) break;//即比一个合数大的质数和该合数的乘积可用一个更大的合数和比其小的质数相乘得到
             }
       }
}

最难理解的是这句话：

if (i % prime[j] == 0) break;
要理解这句话，（顺便不严谨地）证明这个算法的时间复杂度和正确性，要从下面两个方面：

• 每个数至少被访问一次
• 每个数至多被访问一次

• 每个数至少被访问一次

• 对于质数，一定会在i的循环中访问到，并确定为质数。
• 对于质数，一定会在i的循环中访问到，并确定为质数。

• 对于合数，一定可以分解为一个最小素因子和其他数的乘积。

比如

合数 (i = p)（最小素因子）(* a);
若 (i%prime[k] ==0);
则(p * a * prime[k+1]) 可以被后面的 (a * prime[k+1]) 再乘以素数 (p) 筛选出来，
（显而(p<prime[k+1])） 所以(i%prime[k] == 0) 时要停止。
证毕
综上所述，每个数被访问一次且仅访问一次！因此整个算法的复杂度是O(n)的。

$ $

---

$ $

$ $

扩展：米勒罗宾算法

如果我们做题的时候空间不够怎么办？没办法筛素数了怎么办？没事，交给你一个办法!

原理

大于等于5的质数一定和6的倍数相邻

证明：

令(x≥1)，将大于等于(5)的自然数表示如下： ······$6x-1，6x，6x+1，6x+2，6x+3，6x+4，6x+5，6(x+1），6(x+1)+1 ······ $可以看到，不在6的倍数两侧，即(6x)两侧的数为(6x+2，6x+3，6x+4，由于2(3x+1)，3(2x+1)，2(3x+2)，)所以它们一定不是素数，再除去(6x)本身，显然，素数要出现只可能出现在(6x)的相邻两侧。这里要注意的一点是，在(6)的倍数相邻两侧并不是一定就是质数。 根据以上规律，判断质数可以(6)个为单元快进，即将直观判断法的循环中(i++)步长加大为(6)，加快判断速度。

原因是，假如要判定的数为(n)，则(n)必定是(6x-1)或(6x+1)的形式，对于循环中(6i-1，6i，6i+1,6i+2，6i+3，6i+4，)其中如果(n)能被(6i，6i+2，6i+4)整除，则(n)至少得是一个偶数，但是(6x-1)或(6x+1)的形式明显是一个奇数，故不成立；另外，如果(n)能被(6i+3)整除，则(n)至少能被(3)整除，但是(6x)能被(3)整除，故(6x-1)或(6x+1)（即(n)）不可能被(3)整除，故不成立。综上，循环中只需要考虑(6i-1)和(6i+1)的情况，即循环的步长可以定为(6)，每次判断循环变量(k)和(k+2)的情况即可，理论上讲整体速度应该会是直观判断法改进的3倍。

米勒罗宾单次复杂度约为(log(n)*k)((k)为常数且一般取(3)) 判断的素数在(10)亿以内进行(50w)次计算也不会超时

Code：

bool check(int a){
    if(a==1) return 0;
    if(a==2||a==3) return 1;
    if(a%6!=1&&a%6!=5) return 0;
    int temp=sqrt(a);
    for(int i=5;i<=temp;i+=6)
    {
        if(a%i==0||a%(i+2)==0) return 0;
    }
    return 1;
}