入门杜教筛啦。
http://blog.csdn.net/skywalkert/article/details/50500009(好文!)
可以在$O(N^{frac{2}{3}})或O(N^{frac{3}{4}})$的复杂度内解决求某些数论函数f(n)(或f的前缀和S(n)$)的值。
先来看看原理是什么。(接下来推导如何求数论函数f(n)的前缀和S(n))
现在有两个数论函数$f( )和g( )$
(同时定义f的前缀和函数$S(n)=sum_{i=1}^{n}f(i)$)
有狄利克雷乘积可知:
$$f*g(n)=sum_{i|n}f(frac{n}{i})g(i)quad(=sum_{i|n}f(i)g(frac{n}{i}))$$
那么,则有如下结论:
$$sum_{n=1}^{N}f*g(n)=sum_{i=1}^{N}g(i)S(lfloor frac{N}{i} floor)$$
证明如下:
$$egin{align*}
sum_{n=1}^{N}f*g(n)&=sum_{n=1}^{N}sum_{i|n}f(frac{n}{i})g(i)\
&=sum_{i=1}^{N}g(i)sum_{i=1}^{lfloor frac{N}{i}
floor}f(i)\
&=sum_{i=1}^{N}g(i)S(lfloor frac{N}{i}
floor)
end{align*}$$
然后把右边和式里的$g(1)S(N)$那一项提出来得到:
$$g(1)S(N)=sum_{n=1}^{N}f*g(n)-sum_{i=2}^{N}g(i)S(lfloor frac{N}{i} floor)$$
通常令数论函数$g(n)=I(n)=1$(恒等函数$l(n)=1$,完全积性)
到目前为止,上式就是我们进行杜教筛的基础了。
因为左边的S(N)就是答案,而右边同时又可以用分块的方式计算。
不少刚刚入门的同学会有一个疑问,等式右边的后半部分确实可以分块计算,但是前半部分怎么办呢?
其实,一般前面的$sum_{n=1}^{N}f*g(n)$都是可以O(1)计算出来的。
下面来看两个例子:
(一)、求莫比乌斯函数$mu(n)$的前缀和函数$S(n),n leq 10^9$
首先添加一个辅助函数g(n)=l(n)=1,
然后重复上面的过程,可以得到
$$g(1)S(N)=sum_{n=1}^{N}mu*g(n)-sum_{i=2}^{N}g(i)S(lfloor frac{N}{i} floor)$$
$$S(N)=sum_{n=1}^{N}mu*g(n)-sum_{i=2}^{N}S(lfloor frac{N}{i} floor)$$
现在来看看$sum_{n=1}^{N}mu*g(n)$怎么求:
$$egin{aligned}
sum_{n=1}^{N}mu*g(n)&=sum_{n=1}^{N}sum_{i|n}mu(i)g(frac{n}{i})\
&=sum_{n=1}^{N}sum_{i|n}mu(i)\
&=sum_{n=1}^{N}[n==1]\
&=1
end{aligned}$$
上面的化简用到了刚刚学莫比乌斯函数时的一个结论:
$$sum_{i|n}mu(i)=[n==1]$$
到此,我们得到:
$$S(N)=1-sum_{i=2}^{N}S(lfloor frac{N}{i} floor)$$
实现方式是分块计算+记忆化递归处理(用map或者hash表记忆化)
(二). 求欧拉函数$phi(n)$的前缀和函数$S(n),n leq 10^9$
同样地,添加一个辅助函数g(n)=l(n)=1,
然后重复上面的过程,可以得到
$$g(1)S(N)=sum_{n=1}^{N}phi*g(n)-sum_{i=2}^{N}g(i)S(lfloor frac{N}{i} floor)$$
$$S(N)=sum_{n=1}^{N}phi*g(n)-sum_{i=2}^{N}S(lfloor frac{N}{i} floor)$$
$sum_{n=1}^{N}phi*g(n)$又怎样求呢:
$$egin{aligned}
sum_{n=1}^{N}phi*g(n)&=sum_{n=1}^{N}sum_{i|n}phi(i)g(frac{n}{i})\
&=sum_{n=1}^{N}sum_{i|n}phi(i)\
&=sum_{n=1}^{N}n\
&=frac{(1+n)n}{2}
end{aligned}$$
上面的化简用到了这样一个结论:
$$sum_{i|n}phi(i)=n$$
所以我们得到:
$$S(N)=frac{(1+n)n}{2}-sum_{i=2}^{N}S(lfloor frac{N}{i} floor)$$
这个同样是实现方式是分块计算+记忆化递归处理(用map或者hash表记忆化)
下面是代码具体实现:
关于时间复杂度的分析不太会。记了一下结论。
通常不做任何处理,就直接杜教筛的话(分块计算+记忆化递归处理),复杂度是$O(N^{frac{3}{4}})$
但是如果预处理出前$N^{frac{2}{3}}$个前缀和,那么总的复杂度就可以降到$O(N^{frac{2}{3}})$
BZOJ 3944: Sum,杜教筛入门题。
多个询问,给出N,求$sum_{n=1}^{N}mu(n)$和$sum_{n=1}^{N}phi(n)$
也就是求上面的两个例子。
这里直接给出代码,用的是预处理前$N^{frac{2}{3}}$个前缀和+hash表进行记忆化。
复杂度$O(N^{frac{2}{3}})$
#include<cstdio> #include<cstring> #include<iostream> #include<algorithm> #define DJM 1664510 //#define DJM 10 #define ll long long using namespace std; ll phi[DJM+50],mu[DJM+50]; struct Pii{ int x; ll a,b; Pii(int _x=0,ll _a=0,ll _b=0):x(_x),a(_a),b(_b){} }nl; struct Hash_Table{// #define hmod 1425367 int nxt[hmod],head[hmod],hnt; Pii info[hmod]; Hash_Table(){hnt=2;} void Push(Pii rtm){ static int u; u=rtm.x%hmod; info[hnt]=rtm; nxt[hnt]=head[u]; head[u]=hnt++; } Pii Find(int x){ static int u; u=x%hmod; for(int i=head[u];i;i=nxt[i]) if(info[i].x==x) return info[i]; return nl; } }H; void Sieve(){ static bool np[DJM+50]; static int prime[DJM+50],pnt; phi[1]=mu[1]=1; for(int i=2;i<=DJM;i++){ if(!np[i]) prime[++pnt]=i,mu[i]=-1,phi[i]=i-1; for(int j=1;j<=pnt&&i<=DJM/prime[j];j++){ np[i*prime[j]]=1; if(i%prime[j]){mu[i*prime[j]]=-mu[i]; phi[i*prime[j]]=phi[i]*phi[prime[j]];} else{phi[i*prime[j]]=phi[i]*prime[j]; break;} } } for(int i=2;i<=DJM;i++) mu[i]+=mu[i-1],phi[i]+=phi[i-1]; } Pii DJ_Sieve(int x){ if(x<=DJM) return Pii(x,mu[x],phi[x]); if(H.Find(x).x) return H.Find(x); Pii tmp,now=Pii(x,1,(1ll+x)*x/2); for(ll i=2,last;i<=x;i=last+1){ last=x/(x/i); tmp=DJ_Sieve(x/i); now.a-=tmp.a*(last-i+1); now.b-=tmp.b*(last-i+1); } H.Push(now); return now; } int main(){ Sieve(); int Case,n; Pii ans; scanf("%d",&Case); for(int i=1;i<=Case;i++){ scanf("%d",&n); if(n==0) {printf("0 0 "); continue;} ans=DJ_Sieve(n); printf("%lld %lld ",ans.b,ans.a); } return 0; }