杜教筛学习笔记

杜教筛学习笔记

杜教筛被用于求解某一数论函数 \(f\) 的前缀和，即对于形如 \(S(n)=\sum_{i=1}^nf(i)\) 形式的函数 \(S\)，杜教筛能够在小于线性复杂度的复杂度内求解。

算法思想

尝试构造一个函数 \(S\) 的递推式。选择一个数论函数 \(g\)，那么根据狄利克雷卷积可以得到：

\[\begin{aligned}\sum_{i=1}^{n}(f*g)(i)&=\sum_{i=1}^{n}\sum_{d|i}g(d)f(\frac{i}{d})\\&=\sum_{d=1}^{n}g(d)\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}f(i)\\&=\sum_{i=1}^ng(i)S(\lfloor\frac{n}{i}\rfloor)\end{aligned} \]

那么接着我们就可以简单的推得：

\[\begin{aligned}g(1)S(n)&=\sum_{i=1}^ng(i)S(\lfloor\frac{n}{i}\rfloor)-\sum_{i=2}^ng(i)S(\lfloor\frac{n}{i}\rfloor)\\&=\sum_{i=1}^n(f*g)(i)-\sum_{i=2}^ng(i)S(\lfloor\frac{n}{i}\rfloor)\end{aligned} \]

不难看出，只需要能够快速求解 \(\sum_{i=1}^n(f*g)(i)\) 和 \(\sum_{i=2}^ng(i)S(\lfloor\frac{n}{i}\rfloor)\)，就可以快速求出 \(g(1)S(n)\)，也能够快速求得 \(S(n)\)。

事实上，不论 \(f\) 是否是积性函数，只要构造出的函数 \(g\) 能够满足上面的式子的推导，就可以使用杜教筛进行求解。

时间复杂度

在不能预处理的情况下，杜教筛的复杂度是 \(O(n^{\frac{3}{4}})\) 的。

而在可以进行预处理的情况下，进行前 \(n^{\frac{2}{3}}\) 个数的预处理能够让时间复杂度降至 \(O(n^{\frac{2}{3}})\)。

更加值得注意的是，上述时间复杂度全部建立在记忆化的基础上，因此需要采用 map 或 unordered_map 进行存储。

常见积性函数及其卷积形式

\[\varepsilon=\mu*1\\ d=1*1\\ \sigma=\text{id}*1\\ \varphi=\mu*\text{id} \]

其中，\(\varepsilon\) 单位函数，\(1\) 是常函数，\(d\) 是因子函数，\(\sigma\) 是除数函数，\(\text{id}\) 是恒等函数。

代码实现

下面以 \(\mu\) 函数的杜教筛为实例代码。那么所求 \(S(n)=\sum_{i=1}^n\mu(i)\)。

不妨令 \(g=1\)，所以：

\[\begin{aligned}g(1)S(n)&=S(n)\\&=\sum_{i=1}^{n}(\mu*1)(i)-\sum_{i=2}^{n}g(i)S(\lfloor\frac{n}{i}\rfloor)\\&=\sum_{i=1}^n\varepsilon(i)-\sum_{i=2}^n{S(\lfloor\frac{n}{i}\rfloor)}\\&=1-\sum_{i=2}^{n}S(\lfloor\frac{n}{i}\rfloor)\end{aligned} \]

不难看出后半部分可以用数论分块求解。

int cnt;
int p[N+5],mu[N+5];
bool vis[N+5];
unordered_map<int,ll>sum_mu;

void EulerS(){
    mu[1]=1;
    for(int i=2;i<=N;i++){
        if(!vis[i]){
            p[++cnt]=i;
            mu[i]=-1;
        }
        for(int j=1;p[j]<=N/i;j++){
            vis[i*p[j]]=true;
            if(i%p[j]==0){
                mu[i*p[j]]=0;
                break;
            }
            mu[i*p[j]]=-mu[i];
        }
    }
    for(int i=1;i<=N;i++)mu[i]+=mu[i-1];
    return ;
}//线性筛预处理

ll SumMu(ll n){
    if(n<=N)return mu[n];//线性筛预处理出的答案可以不用求解
    if(sum_mu[n])return sum_mu[n];//之前求解过的答案可以不用求解
    ll res=1;
    for(ll l=2,r;l<=n;l=r+1){
        r=n/(n/l);
        res-=SumMu(n/l)*(r-l+1);
    }//求解前缀和
    return sum_mu[n]=res;//记忆化剪枝
}