更快的质数筛: Wheel Factorization

能在短时间内筛到 \(10^9\)。运用了类似分块的思想

如何做

首先,我们取一个合数 \(k\)。先把 \([1,k]\) 之间筛一下。

然后我们顺便得到哪些数和 \(k\) 互质。

有这样一个引理:

对于 \(k\),设 \(Q(i)\) 表示 \(i\) 是否和 \(k\) 互质,是为 \(1\),否为 \(0\)。则:

\(Q(i)=Q(i+k),\forall i\in \N^{*}\)。即 \(Q\) 的周期是 \(k\)

接下来,我们连续枚举长度为 \(k\) 的区间。即,\([k+1,2k],[2k+1,3k],...\)。对于这些区间中的每个数,它们是否和 \(k\) 互质的情况,根据那个引理,我们可以知道。而我们从 \([k+1,2k]\) 开始枚举,它们都 \(>k\)。显然,一个 \(>k\) 的数,如果和 \(k\) 不互质,一定不可能是质数,因此就被筛掉了。

对于每个区间,我们只需要考虑剩下的 \(\varphi(k)\) 个和 \(k\) 互质的就行了。那它们一定没有 \(k\) 中的质因子。也就是说,我们在筛质数的时候,跳过 \(k\) 的质因子,就可以完美遍历到这些和 \(k\) 互质的数了。

因此复杂度是 \(O(n\times \dfrac{\varphi(k)}{k})\)。写一个类似埃氏筛的东西,复杂度就是对的了。

如何让 \(\dfrac{\varphi(k)}{k}\) 尽可能小呢?首先,\(k\) 中的质因子次数超过 \(1\) 肯定是浪费的,无效;然后 \(k\) 的质因子越多越好

对于 \(1e9\) 的范围,取 \(k=2\times 3\times 5\times 7\times 11\times 13\times 17=510510\),事实证明,它就够快了。 此时

\(\dfrac{\varphi(k)}{k}\approx 0.1805\)

模板题:spoj某题

代码

#include <bits/stdc++.h>
using namespace std;
namespace Flandre_Scarlet
{
    #define N 510515
    #define M 181000000
    #define ull unsigned long long
    #define F(i,l,r) for(int i=l;i<=r;++i)
    #define D(i,r,l) for(int i=r;i>=l;--i)
    #define Tra(i,u) for(int i=G.h[u],v=G.to(i);~i;i=G.nx(i),v=G.to(i)) if (i>=0)
    #define MEM(a,x) memset(a,x,sizeof(a))
    #define FK(a) MEM(a,0)
    #define sz(x) ((int)x.size())
    #define all(x) x.begin(),x.end()
    #define p_b push_back
    #define pii pair<int,int>
    #define fir first
    #define sec second
    int K=510510,phiK=92160; // 2*3*5*7*11*13*17, 前7个
    int pr[60000000]; bool notp[N];
    int pos[N],id[N]; // 互质位置, 互质数的编号
    bool tmp[N];
    #define pid(x) ((x/K)*phiK+id[x%K])
    void Init()
    {
        int n=1000000000; // 1e9
        int&c=pr[0],&cp=pos[0];
        notp[1]=1; c=cp=0;
        F(i,2,K) // 先筛一下 [1,K]
        {
            if (!notp[i]) pr[++c]=i;
            for(int j=1;j<=c and pr[j]<=K/i;++j)
            {
                notp[i*pr[j]]=1;
                if (i%pr[j]==0) break;
            }
        }
        F(i,1,K) if (i%2 and i%3 and i%5 and i%7 and i%11 and i%13 and i%17) // 互质
        {
            ++cp;
            pos[cp]=i;
            id[i]=cp;
        }
        F(i,2,1959)
        {
            int r=i*K,l=r-K+1;
            FK(tmp);
            // 每次的tmp[i]表示l+i-1是否不是质数
            // tmp[i]=1表示被筛掉, 最后tmp[i]=0说明是质数
            for(int j=8;j<=c and pr[j]*pr[j]<=r;++j) // pr[j]筛到的最小数是pr[j]*pr[j], 比它小的会更早被筛到
            {
                int fir=pr[j]*max(pr[j],(l-1)/pr[j]+1);
                if (!(fir&1)) fir+=pr[j];
                fir-=(l-1);
                while(fir<=K) tmp[fir]=1,fir+=(pr[j]<<1); // 只需要筛奇数倍, 因为偶数倍含有2, 肯定不互质
            }
            F(j,1,cp) if (!tmp[pos[j]]) pr[++c]=l+pos[j]-1;
        }
        for(int i=1;i<=c and pr[i]<=n;i+=500)
        {
            printf("%d\n",pr[i]);
        }
    }
    void Input()
    {

    }
    void Sakuya()
    {

    }
    void IsMyWife()
    {
        Init();
        Input();
        Sakuya();
    }
}
#undef int //long long
int main()
{
    Flandre_Scarlet::IsMyWife();
    return 0;
}
posted @ 2021-10-12 14:44  Flandre-Zhu  阅读(617)  评论(0编辑  收藏  举报