数学筛法

有的时候怕忘记，写篇小博客记录一下。

线性筛素数

inline void init(int n)
{
    for(int i=2;i<=n;i++)
    {
        if(!vis[i]) prime[++cnt]=i;
        for(int j=1;j<=cnt && i*prime[j];j++)
        {
            vis[i*prime[j]]=1;
            if(!(i%prime[j])) break;
        }
    }
    return;
}

线性筛求欧拉函数

inline void init(int n)
{
    memset(vis,0,sizeof(vis));
    phi[1]=1;
    for(int i=2;i<=n;i++)
    {
        if(!vis[i]) {prime[++cnt]=vis[i]=i，phi[i]=i-1;}
        for(int j=1;j<=cnt && i*prime[j]<=n;j++)
        {
            vis[i*prime[j]]=prime[j];
            phi[i*prime[j]]=phi[i]*(i%prime[j]?prime[j]-1:prime[j]);
        }
    }
    for(int i=1;i<=n;i++) sum[i]=sum[i-1]+phi[i];
    return;
}

线性筛求莫比乌斯函数

inline void init(int n)
{
    mu[1]=1;
    for(int i=2;i<=n;i++)
    {
        if(!vis[i]) {mu[i]=-1;prime[++cnt]=i;}
        for(int j=1;j<=cnt && i*prime[j]<=n;j++)
        {
            vis[i*prime[j]]=1;
            if(!(i%prime[j])) {mu[i*prime[j]]=0;break;}
            else mu[i*prime[j]]=-mu[i];
        }
    }
    for(int i=1;i<=n;i++) sum[i]=sum[i-1]+mu[i];
    return;
}

线性筛求约数个数

用 $d_i$ 表示 $i$ 的约数个数，$num_i$ 表示 $i$ 的最小质因子出现次数。

inline void init(int n)
{
    d[1]=1;
    for(int i=2;i<=n;i++)
    {
        if(!vis[i]) 
        {
            prime[++cnt]=i;
            d[i]=2,num[i]=1;
        }
        for(int j=1;j<=cnt && i*prime[j]<=n;j++)
        {
            vis[i*prime[j]]=1;
            if(!(i%prime[j]))
            {
                num[i*prime[j]]=num[i]+1;
                d[i*prime[j]]=d[i]/num[i*prime[j]]*(num[i*prime[j]]+1);
                break;
            }
            else num[i*prime[j]]=1,d[i*prime[j]]=d[i]*2;
        }
    }
    return;
}

线性筛求约数和

用 $s_i$ 表示 $i$ 的约数和，$g_i$ 表示 $i$ 的最小质因子的 $p^0+p^1+p^2+\dots +p^k$.

inline void init(int n)
{
    s[1]=g[1]=1;
    for(int i=2;i<=n;i++)
    {
        if(!vis[i]) vis[i]=1,prime[++cnt]=i,g[i]=i+1,s[i]=i+1;
        for(int j=1;j<=cnt && i*prime[j]<=n;j++)
        {
            vis[i*prime[j]]=1;
            if(!(i%prime[j]))
            {
                g[i*prime[j]]=g[i]*prime[j]+1;
                s[i*prime[j]]=s[i]/g[i]*g[i*prime[j]];
                break;
            }
            else
            {
                s[i*prime[j]]=s[i]*s[prime[j]];
                g[i*prime[j]]=prime[j]+1;
            }
        }
    }
    return;
}

杜教筛

对于数论函数 $f$，杜教筛可以在低于线性时间的复杂度内计算 $S(n)=\sum\limits_{i=1}^{n}f(i)$。

原理

我们要想办法构造一个 $S(n)$ 关于 $S(\lfloor\frac{n}{i}\rfloor)$ 的递推式。

对于任意一个数论函数 $g$，一定满足：

$$\begin{aligned} \sum\limits_{i=1}^n(f*g)(i)&=\sum\limits_{i=1}^n\sum\limits_{d|i}g(d)f\left(\frac{i}{d}\right)\\ &=\sum\limits_{i=1}^n g(i)S\left(\left\lfloor\frac{n}{i}\right\rfloor\right) \end{aligned}$$

那么可以得到递推式：

$$\begin{aligned} g(1)S(n)&=\sum\limits_{i=1}^n g(i)S\left(\left\lfloor\frac{n}{i}\right\rfloor\right)-\sum\limits_{i=2}^n g(i)S\left(\left\lfloor\frac{n}{i}\right\rfloor\right)\\ &=\sum\limits_{i=1}^n(f*g)(i)-\sum\limits_{i=2}^n g(i)S\left(\left\lfloor\frac{n}{i}\right\rfloor\right) \end{aligned}$$

假如我们可以构造恰当的数论函数 $g$ 使得：

• 可以快速计算 $\sum_{i=1}^n(f * g)(i)$。
• 可以快速计算 g 的单点值，以用数论分块求解 $\sum_{i=2}^ng(i)S\left(\left\lfloor\frac{n}{i}\right\rfloor\right)$。

则我们可以在较短时间内求得 $g(1)S(n)$。

那么就得到 $S(n)=\frac{\sum\limits_{i=1}^n(f*g)(i)-\sum\limits_{i=2}^n g(i)S\left(\left\lfloor\frac{n}{i}\right\rfloor\right)}{g(1)}$。

时间复杂度

反正是 $\Theta(n^{\frac{2}{3}})$

证明

例题1 P4213【模板】杜教筛（Sum）

点击查看代码

#include<bits/stdc++.h>
#define int long long
using namespace std;
const int MAXN=5e6+5;
 
int T,n,cnt,a[MAXN];
int prime[MAXN],vis[MAXN],mu[MAXN],phi[MAXN],summu[MAXN],sumphi[MAXN];
int f[MAXN],tot[MAXN];
map<int,int>ansmu;
map<int,int>ansphi;
 
inline void init(int n)
{
    mu[1]=1; phi[1]=1;
    for(int i=2;i<=n;i++)
    {
        if(!vis[i]) prime[++cnt]=i,mu[i]=-1,phi[i]=i-1;
        for(int j=1;j<=cnt && i*prime[j]<=n;j++)
        {
            vis[i*prime[j]]=1;
            if(!(i%prime[j])) {mu[i*prime[j]]=0;phi[i*prime[j]]=phi[i]*prime[j];break;}
            else mu[i*prime[j]]=-mu[i],phi[i*prime[j]]=phi[i]*(prime[j]-1);
        }
    }
    for(int i=1;i<=n;i++) summu[i]=summu[i-1]+mu[i];
    for(int i=1;i<=n;i++) sumphi[i]=sumphi[i-1]+phi[i];
    return;
}
 
inline int getmu(int x)
{
    if(x<=MAXN) return summu[x];
    if(ansmu[x]) return ansmu[x];
    int ans=1;
    for(int l=2,r;l<=x;l=r+1)
    {
        r=x/(x/l);
        ans-=(r-l+1)*getmu(x/l);
    }
    return ansmu[x]=ans;
}
 
inline int getphi(int x)
{
    if(x<=MAXN) return sumphi[x];
    if(ansphi[x]) return ansphi[x];
    int ans=(x*(x+1))/2ll;
    for(int l=2,r;l<=x;l=r+1)
    {
        r=x/(x/l);
        ans-=(r-l+1)*getphi(x/l);
    }
    return ansphi[x]=ans;
}
 
signed main()
{
    ios_base::sync_with_stdio(false);
    cin.tie(0),cout.tie(0);
    cin>>T; init(MAXN);
    while(T--)
    {
        cin>>n;
        printf("%lld %lld\n",getphi(n),getmu(n));
    }
    return 0;
}