【bzoj3994】 SDOI2015—约数个数和

http://www.lydsy.com/JudgeOnline/problem.php?id=3994 (题目链接)

题意

　　多组询问，给出${n,m}$，求${sum_{i=1}^nsum_{j=1}^m d(i×j)}$，${d(i×j)}$为${ij}$的约数个数。

Solution

　　看到这个式子感觉无从下手，这个${d(ij)}$比较丑，有一个比较经典的公式：$${d(nm)=sum_{i|n}sum_{j|m} [gcd(i,j)=1]}$$

　　这个是怎么得来的呢。每一个${nm}$的约数都可以表示成这样的形式：${i×frac{m}{j}}$。其中${i}$是${n}$的约数，${j}$是${m}$的约数。

　　那么如果我们直接枚举${i,j}$来统计的话，可能会有重复，所以就要加上一个条件：${[gcd(i,j)=1]}$，也就是${i,j}$互质。那他们不互质为什么就不行呢，假设${i,j}$都含有一个约数${p}$，那么${frac{m}{j}}$就表示从${m}$中拿掉${j}$用剩下的去组成约数，也就是从约数中拿掉了${p}$；而${i}$就是从${n}$中拿${i}$去组成约数，也就向约数中加入了${p}$。那么我们拿掉一个${p}$又加入一个${p}$，不是吃多了吗→_→，这样自然会算重复啦。

　　解决了这个问题，我们就很好做了。开始推式子。

egin{aligned}  sum_{i=1}^{n}sum_{j=1}^{m}d(ij) =& sum_{i=1}^nsum_{j=1}^msum_{u|i}sum_{v|j}[gcd(u,v)=1]  \  =&sum_{u=1}^nsum_{v=1}^m[gcd(u,v)=1]lfloorfrac{n}{i} floorlfloorfrac{m}{j} floor  \  =&sum_{t=1}^nμ(t)sum_{i=1}^{lfloor{n/t} floor}sum_{j=1}^{lfloor{m/t} floor}lfloorfrac{n}{ti} floorlfloorfrac{m}{tj} floor  \  =&sum_{t=1}^nμ(t)sum_{i=1}^{lfloor{n/t} floor}lfloorfrac{lfloor{n/t} floor}{i} floorsum_{j=1}^{lfloor{m/t} floor}lfloorfrac{lfloor{m/t} floor}{j} floor   end{aligned}

　　我们令${f(n)=sum_{i=1}^nlfloor{n/i} floor}$，那么${f}$是可以分段${O(nsqrt{n})}$的预处理出来的。$${原式=sum_{t=1}^{n}μ(t)f(lfloor{n/t} floor)f(lfloor{m/t} floor)}$$

　　这样我们就可以预处理出${μ}$的前缀和，然后分段求解答案就可以了。

细节

　　LL，预处理不要太多会TLE。

代码

// bzoj3994
#include<algorithm>
#include<iostream>
#include<cstdlib>
#include<cstring>
#include<cstdio>
#include<vector>
#include<cmath>
#define LL long long
#define inf 2147483640
#define Pi acos(-1.0)
#define free(a) freopen(a".in","r",stdin),freopen(a".out","w",stdout);
using namespace std;

const int maxn=50010;
LL n,m,s[maxn],mu[maxn],f[maxn];
int p[maxn],vis[maxn];

int main() {
	int T;scanf("%d",&T);
	s[1]=mu[1]=1;
	for (int i=2;i<maxn;i++) {
		if (!vis[i]) p[++p[0]]=i,mu[i]=-1;
		for (int j=1;j<=p[0] && p[j]*i<maxn;j++) {
			vis[i*p[j]]=1;
			if (i%p[j]==0) {mu[i*p[j]]=0;break;}
			else mu[i*p[j]]=-mu[i];
		}
		s[i]=s[i-1]+mu[i];
	}
	for (int i=1;i<maxn;i++)
		for (int j=1,k;j<=i;j=k+1) {
			k=i/(i/j);
			f[i]+=(LL)(k-j+1)*(i/j);
		}
	while (T--) {
		scanf("%lld%lld",&n,&m);
		if (n>m) swap(n,m);
		LL ans=0;
		for (int i=1,j;i<=n;i=j+1) {
			j=min(n/(n/i),m/(m/i));
			ans+=f[n/i]*f[m/i]*(s[j]-s[i-1]);
		}
		printf("%lld
",ans);
	}
	return 0;
}