洛谷P3327 [SDOI2015]约数个数和 【莫比乌斯反演】

题目

设d(x)为x的约数个数,给定N、M,求(sum_{i = 1}^{N} sum_{j = 1}^{M} d(ij))

输入格式

输入文件包含多组测试数据。第一行,一个整数T,表示测试数据的组数。接下来的T行,每行两个整数N、M。

输出格式

T行,每行一个整数,表示你所求的答案。

输入样例

2
7 4
5 6

输出样例

110
121

提示

1<=N, M<=50000

1<=T<=50000

题解

好神的题【是我太弱吧
首先上来就伤结论。。

题目所求

[ans = sum_{i = 1}^{N} sum_{j = 1}^{M} d(ij) ]

有一个这样的结论:

[d(ij) = sum_{x|i}sum_{y|j} [gcd(x,y) == 1] ]

那么就转化为了:

[ans =sum_{i = 1}^{N} sum_{j = 1}^{M} sum_{x|i}sum_{y|j} [gcd(x,y) == 1] ]

我们考虑对于每一对互质的x、y,x会被枚举(lfloor frac{N}{x} floor)次,y会被枚举(lfloor frac{M}{y} floor)
所以有

[ans =sum_{i = 1}^{N} sum_{j = 1}^{M} lfloor frac{N}{i} floor lfloor frac{M}{j} floor [gcd(i,j) == 1] ]

那么可以进行莫比乌斯反演了

[f(n) = sum_{i = 1}^{N} sum_{j = 1}^{M} lfloor frac{N}{i} floor lfloor frac{M}{j} floor [gcd(i,j) == n] ]

[F(n) = sum_{i = 1}^{N} sum_{j = 1}^{M} lfloor frac{N}{i} floor lfloor frac{M}{j} floor [n | gcd(i,j)] ]

那么有

[egin{aligned} F(d) &= (sum_{i = 1}^{N}lfloor frac{N}{i} floor) * (sum_{j = 1}^{M} lfloor frac{M}{j} floor) [d | gcd(i,j)] \ &= (sum_{i = 1}^{lfloor frac{N}{d} floor}lfloor frac{N}{id} floor) * (sum_{j = 1}^{lfloor frac{M}{d} floor} lfloor frac{M}{jd} floor) \ &= (sum_{i = 1}^{lfloor frac{N}{d} floor}lfloor frac{lfloor frac{N}{d} floor}{i} floor) * (sum_{j = 1}^{lfloor frac{M}{d} floor} lfloor frac{lfloor frac{M}{d} floor}{j} floor) end{aligned} ]

其中(sum_{i = 1}^{N}lfloor frac{N}{i} floor)可以(O(nsqrt{n}))预处理出,我们记为(sum(n))
那么

[F(n) = sum(lfloor frac{N}{n} floor) * sum(lfloor frac{M}{n} floor) ]

[ans = f(1) = sum_{d = 1}^{N} mu(d) * F(d) = sum_{d} mu(d) sum(lfloor frac{N}{d} floor) * sum(lfloor frac{M}{d} floor) ]

分块计算
复杂度(O(Tsqrt{N} + Nsqrt{N}))

#include<iostream>
#include<cmath>
#include<cstdio>
#include<cstring>
#include<algorithm>
#define LL long long int
#define REP(i,n) for (int i = 1; i <= (n); i++)
#define Redge(u) for (int k = h[u],to; k; k = ed[k].nxt)
#define BUG(s,n) for (int i = 1; i <= (n); i++) cout<<s[i]<<' '; puts("");
using namespace std;
const int maxn = 50005,maxm = 100005,INF = 1000000000;
inline int read(){
	int out = 0,flag = 1; char c = getchar();
	while (c < 48 || c > 57) {if (c == '-') flag = -1; c = getchar();}
	while (c >= 48 && c <= 57) {out = (out << 3) + (out << 1) + c - '0'; c = getchar();}
	return out * flag;
}
int prime[maxn],primei,mu[maxn],isn[maxn];
LL sum[maxn];
void init(){
	mu[1] = 1;
	for (int i = 2; i < maxn; i++){
		if (!isn[i]) prime[++primei] = i,mu[i] = -1;
		for (int j = 1; j <= primei && i * prime[j] < maxn; j++){
			isn[i * prime[j]] = true;
			if (i % prime[j] == 0){mu[i * prime[j]] = 0; break;}
			mu[i * prime[j]] = -mu[i];
		}
	}
	for (int i = 1; i < maxn; i++) mu[i] += mu[i - 1];
	for (int n = 1,nxt; n <= 50000; n++){
		for (int i = 1; i <= n; i = nxt + 1){
			nxt = n / (n / i);
			sum[n] += (LL)(nxt - i + 1) * (n / i);
		}
	}
}
int main(){
	init();
	int T = read(),n,m;
	while (T--){
		n = read(); m = read();
		if (n > m) swap(n,m);
		LL ans = 0; int nxt;
		for (int i = 1; i <= n; i = nxt + 1){
			nxt = min(n / (n / i),m / (m / i));
			ans += sum[n / i] * sum[m / i] * (mu[nxt] - mu[i - 1]);
		}
		printf("%lld
",ans);
	}
	return 0;
}

原文地址:https://www.cnblogs.com/Mychael/p/8361823.html