数学筛法

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

线性筛素数

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;
}
posted @ 2023-09-27 16:15  Code_AC  阅读(6)  评论(0编辑  收藏  举报