杜教筛学习笔记

杜教筛可以干什么？他可以在低于 \(O(n)\) 的时间内求出一个积性函数的前缀和。

首先我们要先知道狄利克雷卷积是啥，函数 \(f(n)\) 和函数 \(g(n)\) 的狄利克雷卷积记作 \((f * g)(n)\)，其值为 \((f * g)(n) = \sum_{d \mid n} f(d) g(\frac{n}{d})\)。有很多奇怪的函数与函数之间的狄利克雷卷积是另一个奇怪的函数，比如 \((\mu * 1) = \epsilon,(\varphi * 1) = \textrm{id}\)。

接下来进入杜教筛。首先我们先列出我们要求的函数 \(f(n)\) 和另一个积性函数 \(g(n)\) 的狄利克雷卷积的前缀和：

\[\sum_{i=1}^{n} (f * g)(i) \]

\[=\sum_{i=1}^{n}\sum_{d \mid i}f(d)g(\frac{i}{d}) \]

\[=\sum_{d=1}^{n}g(d)\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}f(i) \]

可以发现这个式子后面的东西是前缀和，记 \(f\) 的前缀和为 \(S\)。

\[=\sum_{d=1}^{n}g(d)S(\lfloor\frac{n}{d}\rfloor) \]

我们尝试把 \(S(n)\) 从右边挪到左边。

这时候我们相减可以得到：

\[g(1)S(n)=\sum_{d=1}^{n}g(d)S(\lfloor\frac{n}{d}\rfloor)-\sum_{d=2}^{n}g(d)S(\lfloor\frac{n}{d}\rfloor) \]

\[=\sum_{i=1}^{n} (f * g)(i)-\sum_{d=2}^{n}g(d)S(\lfloor\frac{n}{d}\rfloor) \]

接下来就可以得出：

\[S(n)=\frac{\sum_{i=1}^{n} (f * g)(i)-\sum_{d=2}^{n}g(d)S(\lfloor\frac{n}{d}\rfloor)}{g(1)} \]

只要我们能制定一个 \(g\) 使得可以快速的求出 \((f * g)\) 的前缀和与 \(g\) 的前缀和即可，上式中分子后半部分可以整除分块。注意我们可以线性筛先筛出较小的 \(n\) 的答案，防止后续递归过慢，这个 \(n\) 可以取 \(N^{(\frac{2}{3})}\)。

这样就搞定了，接下来的目标就是制定一个合理的 \(g\)。

我们给出两个模板题中的例子和一个其他例子。

第一种，若 \(f\) 为莫比乌斯函数。

我们在上文介绍狄利克雷卷积时说过，\((\mu * 1) = \epsilon\)，其中 \(\epsilon\) 仅在 \(n=1\) 时为 \(1\)。这时候 \(g,(f * g)\) 分别为 \(1,\epsilon\)，他们俩的前缀和很好求。

第二种，若 \(f\) 为欧拉函数。

上文也提到过，\((\varphi * 1) = \textrm{id}\)，其中 \(\textrm{id}(n)\) 即为 \(n\)。这时候 \(g,(f * g)\) 分别为 \(1,\textrm{id}\)，前缀和依然好求。

第三种，若 \(f(x)=\varphi(x)\times x^2\)。

乘乘倍了一下，我们考虑 \(g = \textrm{id}^2\)，此时 \((f * g) = \sum_{d \mid n}f(d)g(\frac{n}{d})=\sum_{d \mid n}\varphi(i)\times n^2 = n^3\)。

现在我们的 \(g,(f * g)\) 分别为 \(\textrm{id}^2,\textrm{id}^3\)，前缀和也是好求的。

这玩意时间复杂度不好证，但是在线性筛预处理的 \(n\) 为 \(N^{(\frac{2}{3})}\) 时时间复杂度最优为 \(O(N^{(\frac{2}{3})})\)。

下面给出上面两个例子，也就是 P4213 【模板】杜教筛 的代码。

#include<iostream>
#include<cstdio>
#include<unordered_map>
#define int long long
using namespace std;
const int N=10000000;
int t,n,p[10000005],k,tot,u[10000005],phi[10000005],sumu[10000005],sumphi[10000005];
bool isp[10000005];
void getprime(){
    for(int i=2;i<=10000000;i++) isp[i]=1;u[1]=1;phi[1]=1;
    for(int i=2;i<=10000000;i++){
        if(isp[i]) p[++tot]=i,u[i]=-1,phi[i]=i-1;
        for(int j=1;j<=tot&&i*p[j]<=10000000;j++){
            isp[i*p[j]]=0;
            u[i*p[j]]=u[i]*u[p[j]];
            phi[i*p[j]]=phi[i]*phi[p[j]];
            if(i%p[j]==0){
                u[i*p[j]]=0;
                phi[i*p[j]]=phi[i]*p[j];
                break;
            }
        }
    }
    for(int i=1;i<=10000000;i++){
        sumu[i]=sumu[i-1]+u[i];
        sumphi[i]=sumphi[i-1]+phi[i];
    }
}
unordered_map<int,int>mpu,mpphi;
int getu(int x){
    if(x<=10000000) return sumu[x];
    if(mpu.count(x)) return mpu[x];
    int pre=1,l=2,r=0;
    while(l<=x){
        r=x/(x/l);
        pre-=(r-l+1)*getu(x/l);
        l=r+1;
    }
    mpu[x]=pre;
    return pre;
}
int getphi(int x){
    if(x<=10000000) return sumphi[x];
    if(mpphi.count(x)) return mpphi[x];
    int pre=x*(x+1)/2,l=2,r=0;
    while(l<=x){
        r=x/(x/l);
        pre-=(r-l+1)*getphi(x/l);
        l=r+1;
    }
    mpphi[x]=pre;
    return pre;
}
signed main(){
    cin>>t;
    getprime();
    while(t--){
        cin>>n;cout<<getphi(n)<<' '<<getu(n)<<endl;
    }
    return 0;
}