●BZOJ 3994 [SDOI2015]约数个数和
题链:
http://www.lydsy.com/JudgeOnline/problem.php?id=3994
题解:
莫比乌斯反演
(先定义这样一个符号[x],如果x为true,则[x]=1,否则[x]=0)
首先有这么一个结论:
令d(x)表示x的约数的个数,那么
$d(nm)=\sum_{i|n}\sum_{j|m}[gcd(i,j)==1]$
证明:
设$n=p1^{x1}p2^{x2}p3^{x3}\cdots pk^{xk},m=p1^{y1}p2^{y2}p3^{y3}\cdots pk^{yk}$
则$nm=p1^{x1+y1}p2^{x2+y2}p3^{x3+y3}\cdots pk^{xk+yk}$
由约数定理,$d(nm)=(x1+y1+1)(x2+y2+1)(x3+y3+1)\cdots(xk+yk+1)$
再设$i=p1^{a1}p2^{a2}p3^{a3}\cdots pk^{ak},j=p1^{b1}p2^{b2}p3^{b3}\cdots pk^{bk}$
如果gcd(i,j)=1,那么必须满足a1==0或者b1==0,
如果a1==0,则b1有y1种取值,如果b1==0,则a1有x1种取值,同时a1和b1还可以同时为0
那么就有(x1+y1+1)种情况,
即只考虑p1的指数,就有(x1+y1+1)种情况,同时枚举的i,j如果互质的话,后面的a2,b2,a3,b3...也满足这个条件,
所以满足条件的i,j的对数为$\prod (x_i+y_i+1)$ 和约数定理的形式相同。
有了这个结论,我们来化一化求ANS的式子
$ANS=\sum_{n=1}^{N}\sum_{m=1}^{M}d(nm)$
$\quad\quad=\sum_{n=1}^{N}\sum_{m=1}^{M}\sum_{i|n}\sum_{j|m}[gcd(i,j)==1]$
$\quad\quad=\sum_{i=1}^{N}\sum_{j=1}^{M}[gcd(i,j)==1]\lfloor \frac{N}{i} \rfloor \lfloor \frac{M}{j} \rfloor$
同时由于刚刚入门mobius时,有这么一个式子:
$w(x)=\sum_{d|x}\mu(d)$ 若x==1则w(x)=1,否则w(x)=0
所以:$[gcd(i,j)==1]=\sum_{d|gcd(i,j)}\mu(d)$
那么继续:
$ANS=\sum_{i=1}^{N}\sum_{j=1}^{M}\lfloor \frac{N}{i} \rfloor \lfloor \frac{M}{j} \rfloor\sum_{d|gcd(i,j)}\mu(d)$
$\quad\quad=\sum_{d=1}^{min(n,m)}\mu(d)\sum_{i=1}^{N/d}\lfloor \frac{N}{id}\rfloor\sum_{j=1}^{M/d}\lfloor \frac{M}{jd}\rfloor$
令$f(x)=\sum_{i=1}^{x}\lfloor \frac{x}{i}\rfloor$
则$ANS=\sum_{d=1}^{min(n,m)}\mu(d)f(\lfloor \frac{N}{d} \rfloor)f(\lfloor \frac{M}{d} \rfloor)$
而f(x)就是最开始的d(x)的前缀和。。。但是需要预处理的x的范围小了很多,可以用线筛完成。
代码:
#include<cstdio> #include<cstring> #include<iostream> #define ll long long #define MAXN 50050 using namespace std; ll f[MAXN]; int mu[MAXN]; void Sieve(){ static bool np[MAXN]; static int prime[MAXN],pnt; f[1]=mu[1]=1; for(int i=2,tmp,d;i<=50000;i++){ if(!np[i]) prime[++pnt]=i,mu[i]=-1,f[i]=2; for(int j=1;j<=pnt&&i<=50000/prime[j];j++){ np[i*prime[j]]=1; tmp=i; d=1; while(tmp%prime[j]==0) tmp/=prime[j],d++; f[i*prime[j]]=f[tmp]*(d+1); if(i%prime[j]) mu[i*prime[j]]=-mu[i]; else break; } } for(int i=2;i<=50000;i++) f[i]+=f[i-1],mu[i]+=mu[i-1]; } int main(){ Sieve(); ll ans; int Case,n,m,mini; scanf("%d",&Case); while(Case--){ scanf("%d%d",&n,&m); mini=min(n,m); ans=0; for(int d=1,last;d<=mini;d=last+1){ last=min(n/(n/d),m/(m/d)); ans+=(mu[last]-mu[d-1])*f[n/d]*f[m/d]; } printf("%lld\n",ans); } return 0; }
Do not go gentle into that good night.
Rage, rage against the dying of the light.
————Dylan Thomas