标签：int VID mmu 关于 www. 预处理 stdout git mes
杜教筛真的是一个十分著名的筛法，它那玄学的\(O(n^{\frac23})\)时间复杂度真的是十分神奇。
它主要用途是求积性函数的前缀和（当然，根据差分思想，求一段区间内的值之和也是很简单的）。
假设我们要求的是\(\sum_{i=1}^nf(i)\)（注意，\(f\)是一个积性函数）。
则第一步是寻找（或构造）两个函数\(h\)和\(g\)使得\(h=f*g\)（通常我都找不到两个这样的函数），其中\(*\)表示狄利克雷卷积。
然后我们可以进行一波分析：
\[\sum_{i=1}^nh(i)=\sum_{i=1}^n\sum_{d|i}g(d)·f(\frac id)\]
转化一下就是这个样子：
\[\sum_{i=1}^nh(i)=\sum_{d=1}^ng(d)·\sum_{i=1}^{\lfloor\frac nd\rfloor}f(i)\]
如果用\(sum(n)\)来表示\(\sum_{i=1}^nf(i)\)，那么就可以将式子转化成这个样子：
\[\sum_{i=1}^nh(i)=\sum_{d=1}^ng(d)·sum(\lfloor\frac nd\rfloor)\]
如果我们将右边式子中\(d=1\)的情况单独提出，就可以得到这样一个式子：
\[\sum_{i=1}^nh(i)=\sum_{d=2}^ng(d)·sum(\lfloor\frac nd\rfloor)+g(1)·sum(n)\]
移项可得：
\[g(1)sum(n)=\sum_{i=1}^nh(i)-\sum_{d=2}^ng(d)·sum(\lfloor\frac nd\rfloor)\]
两边同时除以\(g(1)\)可以得到：
\[sum(n)=\frac{\sum_{i=1}^nh(i)-\sum_{d=2}^ng(d)·sum(\lfloor\frac nd\rfloor)}{g(1)}\]
而\(sum(n)\)就是我们所要求的。
所以就不难发现，在\(h\)函数能够快速求出时（如\(e,I,id\)等函数就是非常好的选择），我们就能用杜教筛，在\(O(n^{\frac23})\)的时间内求出\(sum(n)\)。
杜教筛采取了递归的实现方式，而且还使用记忆化加以优化。
根据上面那个式子，我们可以先预处理出\(g\)函数的前缀和，然后对\(sum(\lfloor\frac nd\rfloor)\)进行除法分块。
而除法分块的过程中，我们需要再一次调用\(sum\)函数，这就是一个递归的过程。
在\(n\)比较小的时候，我们可以用线性筛预处理出答案。
于是，就可以通过这种方式来实现我们复杂度为\(O(n^{\frac23})\)的杜教筛了。
注意，记忆化的过程可以使用\(map\)，当然也能使用传说中无比神奇\(O(1)\)查询的\(unordered\_map\)，但是据\(hl666\)大佬的亲身实践，貌似会跑得更慢... ...
这里有一道模板题：【洛谷4213】【模板】杜教筛（Sum）。
大致题意就是让我们求\(\sum_{i=1}^n\phi(i)\)和\(\sum_{i=1}^n\mu(i)\)，而且\(N\le2^{31}-1\)，不能用线性筛。
所以我们要用杜教筛。
首先考虑如何筛\(\sum_{i=1}^n\phi(i)\)。
关于\(\phi\)函数，有这样一个式子\(\phi*I=id\)。
而\(id\)函数是很好求的（\(id(n)=n\)）。
于是可以得到：
\[\sum_{i=1}^n\phi(i)=\frac{\sum_{i=1}^nid(i)-\sum_{i=2}^nI(d)·(\sum_{i-1}^{\lfloor\frac nd\rfloor}\phi(i))}{I(1)}\]
将\(I\)和\(id\)两个函数的值代入，就可以得到一个无比简单的式子：
\[\sum_{i=1}^n\phi(i)=\frac{n(n+1)}2-\sum_{d=2}^n\sum_{i=1}^{\lfloor\frac nd\rfloor}\phi(i)\]
筛\(\sum_{i=1}^n\mu(i)\)也是同理的。
我们可以使用式子\(\mu*I=e\)，化简得到：
\[\sum_{i=1}^n\mu(i)=1-\sum_{d=2}^n\sum_{i=1}^{\lfloor\frac nd\rfloor}\mu(i)\]
这样转化之后就可以用杜教筛来筛了。
代码如下：
#include<bits/stdc++.h>
#define ten(x) (((x)<<3)+((x)<<1))
#define LL long long
using namespace std;
int n;
class FIO
{
private:
#define Fsize 100000
#define tc() (A==B&&(B=(A=Fin)+fread(Fin,1,100000,stdin),A==B)?EOF:*A++)
#define pc(ch) (FoutSize<Fsize?Fout[FoutSize++]=ch:(fwrite(Fout,1,Fsize,stdout),Fout[(FoutSize=0)++]=ch))
int Top,FoutSize;char ch,*A,*B,Fin[Fsize],Fout[Fsize],Stack[Fsize];
public:
inline void read(int &x) {x=0;while(!isdigit(ch=tc()));while(x=ten(x)+(ch&15),isdigit(ch=tc()));}
inline void write(LL x) {if(!x) return (void)(pc('0'));if(x<0) pc('-'),x=-x;while(x) Stack[++Top]=x%10+48,x/=10;while(Top) pc(Stack[Top--]);}
inline void write_char(char x) {pc(x);}
inline void clear() {fwrite(Fout,1,FoutSize,stdout);}
}F;
class Class_DuSieve//杜教筛
{
private:
#define Size 5000000
int Prime_cnt,Prime[Size+5],phi[Size+5],mu[Size+5],SumMu[Size+5];LL SumPhi[Size+5];bool IsNotPrime[Size+5];
map<LL,LL> MemmoryPhi;map<int,int> MemmoryMu;//用map存储答案
public:
Class_DuSieve()//预处理
{
register int i,j;
for(phi[1]=mu[1]=1,i=2;i<=Size;++i)//先用一遍线性筛预处理
{
if(!IsNotPrime[i]) Prime[++Prime_cnt]=i,phi[i]=i-1,mu[i]=-1;
for(j=1;j<=Prime_cnt&&1LL*i*Prime[j]<=Size;++j)
if(IsNotPrime[i*Prime[j]]=1,i%Prime[j]) phi[i*Prime[j]]=phi[i]*(Prime[j]-1),mu[i*Prime[j]]=-mu[i];else {phi[i*Prime[j]]=phi[i]*Prime[j];break;}
}
for(i=1;i<=Size;++i) SumPhi[i]=SumPhi[i-1]+phi[i],SumMu[i]=SumMu[i-1]+mu[i];//求前缀和
}
inline LL GetPhi(int x)//求phi的前缀和
{
if(x<=Size) return SumPhi[x];//如果已经线性筛求出过答案，就返回已经求出的答案
if(MemmoryPhi[x]) return MemmoryPhi[x];//如果记忆化过，返回答案
register int l,r;register LL ans=1LL*x*(x+1)>>1;
for(l=2;l<=x;l=r+1) r=x/(x/l),ans-=(r-l+1)*GetPhi(x/l);//除法分块，递归调用
return MemmoryPhi[x]=ans;//记忆化
}
inline int GetMu(int x)//求mu的前缀和（和上面的函数基本上一样）
{
if(x<=Size) return SumMu[x];
if(MemmoryMu[x]) return MemmoryMu[x];
register int l,r,ans=1;
for(l=2;l<=x;l=r+1) r=x/(x/l),ans-=(r-l+1)*GetMu(x/l);
return MemmoryMu[x]=ans;
}
}DuSieve;
int main()
{
register int T;F.read(T);
while(T--) F.read(n),F.write(DuSieve.GetPhi(n)),F.write_char(' '),F.write(DuSieve.GetMu(n)),F.write_char('\n');
return F.clear(),0;
}标签：int VID mmu 关于 www. 预处理 stdout git mes
原文地址：https://www.cnblogs.com/chenxiaoran666/p/DuSieve.html
