杜教筛
设现在要求积性函数 \(f\) 的前缀和, 设 \(\sum_{i=1}^{n}f(i)=S(n)\)
找一个积性函数\(g\),考虑它们的狄利克雷卷积的前缀和
先枚举 \(d\) 提出 \(g\)
\[\sum_{i=1}^{n} (f*g)(i)= \sum_{d=1}^{n}g(d)\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}f(i) \]把\(\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\)替换为\(S(\lfloor\frac{n}{d}\rfloor)\)
\[\sum_{i=1}^{n} (f*g)(i)= \sum_{d=1}^{n}g(d)S(\lfloor\frac{n}{d}\rfloor) \]可以发现从\(1\)开始的前缀和减去从\(2\)开始的前缀和就是第一项,即
\[g(1)S(n)=\sum_{i=1}^{n}g(i)S(\lfloor\frac{n}{i}\rfloor)-\sum_{i=2}^{n}g(i)S(\lfloor\frac{n}{i}\rfloor) \]所以得到杜教筛的核心式子:
\[g(1)S(n)=\sum_{i=1}^{n}(f*g)(i)-\sum_{i=2}^{n}g(i)S(\lfloor\frac{n}{i}\rfloor) \]则如果可以找到一个合适的积性函数\(g\),使得可以快速算出\(\sum_{i=1}^{n}(f*g)(i)\)和\(g\)的前缀和,便可以用数论分块递归地求解
代码框架如下
il int GetSum_f(cs int n){ // 算 f 前缀和的函数
ri int as=sum_f_g(n); // sum_f_g(n)是算 f * g 的前缀和
// 数论分块
for(ri int l=2,r=0;l<=n;l=r+1){ // 注意从 2 开始
r=(n/(n/l)); // g_sum 是 g 的前缀和,递归求解
as-=(sum_g(r)-sum_g(l-1))*GetSum_f(n/l);
}
return as;
}
复杂度\(O(n^{\frac{3}{4}})\)
还可以进一步优化,即先线性筛出前 \(m\) 个答案,之后再用杜教筛,复杂度是\(O(\frac{n}{\sqrt{n}})\)
当\(m=n^{\frac{2}{3}}\)时,复杂度为\(O(n^\frac{2}{3})\)
可以使用哈希表来存下已经求过的答案,第一篇题解说可以不用。
考虑到上面的求和过程中出现的都是 \(\lfloor\frac{n}{i}\rfloor\)。开一个大小为两倍 \(\sqrt{n}\) 的数组 \(f\) 记录答案。如果现在需要求出 GetSum_f(x) ,若 \(x\le\sqrt{n}\),返回 \(f[x]\) ,否则返回 \(f[sqrt(n)+n/i]\) 即可。
code
求\(\sum_{i=1}^n\varphi(i)\),\(\sum_{i=1}^n\mu(i)\)
题目链接
题解
AC·code
#include<bits/stdc++.h>
#define il inline
#define cs const
#define ri register
#define ll long long
using namespace std;
namespace Q{
il int rd(){
ri int x=0,f=1;ri char c=getchar();
while(c<'0'||c>'9') f=(c=='-')?-1:1,c=getchar();
while(c>='0'&&c<='9') x=x*10+(c^48),c=getchar();
return x*f;
}
il void wt(int x){
if(x<0) x=-x,putchar('-');
if(x>=10) wt(x/10);
return putchar(x%10+48),void();
}
il void wl(ll x){
if(x<0) x=-x,putchar('-');
if(x>=10) wl(x/10);
return putchar(x%10+48),void();
}
}
namespace djs{
cs int N=5e6;
int ss[N],mu[N+5],cnt;
ll fi[N+5];
il void init(){
mu[1]=fi[1]=1;
for(ri int i=2;i<=N;++i){
if(!fi[i]) ss[++cnt]=i,fi[i]=i-1,mu[i]=-1;
for(ri int j=1;j<=cnt&&1ll*ss[j]*i<=N;++j){
if(i%ss[j]){
fi[i*ss[j]]=fi[i]*fi[ss[j]];
mu[i*ss[j]]=-mu[i];
}
else{
fi[i*ss[j]]=fi[i]*ss[j];
break;
}
}
}
for(ri int i=2;i<=N;++i){
fi[i]+=fi[i-1],mu[i]+=mu[i-1];
}
return;
}
unordered_map<int,int> Mu;
unordered_map<int,ll> Fi;
il ll varphi(int n){
if(n<=N) return fi[n];
else if(Fi.count(n)) return Fi[n];
ri ll phi=(1ll+n)*n/2;
for(ri unsigned int l=2,r=0;l<=n;l=r+1){
r=n/(n/l),phi-=1ll*(r-l+1)*varphi(n/l);
}
return Fi[n]=phi;
}
il int get_mu(int n){
if(n<=N) return mu[n];
else if(Mu.count(n)) return Mu[n];
ri int smu=1;
for(ri unsigned int l=2,r=0;l<=n;l=r+1){
r=n/(n/l),smu-=(r-l+1)*get_mu(n/l);
}
return Mu[n]=smu;
}
} using namespace djs;
signed main(){
init();
for(ri int i=Q::rd(),n;i;--i){
n=Q::rd();
Q::wl(varphi(n)),putchar(' ');
Q::wt(get_mu(n)),putchar('\n');
}
return 0;
}