bzoj3994/ SDOI2015 约数个数和（莫比乌斯函数

最新推荐文章于 2021-03-21 18:31:13 发布

原创最新推荐文章于 2021-03-21 18:31:13 发布 · 282 阅读

1 ·

CC 4.0 BY-SA版权

莫比乌斯反演专栏收录该内容

6 篇文章

订阅专栏

本文介绍了一种利用数论和算法技巧解决特定数学问题的方法，即求解给定N和M时，所有i从1到N，j从1到M的i*j的约数个数之和。通过引入莫比乌斯函数和分块思想，将复杂度降低到O(n√n + T√n)。

题目描述

设d(x)为x的约数个数，给定N、M，求 $∑i=1n∑j=1md(ij)\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}d(ij)$

输入格式

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

输出格式

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

输入输出样例

输入 #1

2
7 4
5 6

输出 #1

110
121
说明/提示
1<=N, M<=50000

1<=T<=50000

分析

首先有一个结论（~~别看我我不会证~~
$=\sum\limits_{x|i}\sum\limits_{y|j}[gcd(x,y)==1]$
于是
$\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}d(ij)= \sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}\sum\limits_{x|i}\sum\limits_{y|j}[gcd(x,y)==1]=\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}\sum\limits_{x|i}\sum\limits_{y|j}\sum\limits_{k|x,k|y}\mu(k)=\sum\limits_{k=1}^{n}\mu(k)\sum\limits_{i=1}^{{\tiny\left\lfloor\dfrac{n}{k}\right\rfloor}}\left\lfloor\dfrac{n}{kd}\right\rfloor\sum\limits_{j=1}^{{\tiny\left\lfloor\dfrac{m}{k}\right\rfloor}}\left\lfloor\dfrac{m}{kj}\right\rfloor$
记 $\sum\limits_{i=1}^{n}\left\lfloor\dfrac{n}{i}\right\rfloor$ ，于是
$\sum\limits_{k=1}^{n}\mu(k)*g(\left\lfloor\dfrac{n}{k}\right\rfloor)*g(\left\lfloor\dfrac{m}{k}\right\rfloor)$ ，除法分块即可。
总复杂度 $O(nn+Tn)O(n\sqrt{n}+T\sqrt{n})$

代码如下

#include <bits/stdc++.h>
#define N 50005
#define LL long long
using namespace std;
int g[N], mu[N], s[N], x[N], p[N], cnt, maxn = 50000;
LL ans, z = 1;
int main(){
	int i, j, n, m, t, l, r, T;
	mu[1] = 1;
	for(i = 2; i <= maxn; i++){
		if(!x[i]) x[i] = 1, p[++cnt] = i, mu[i] = -1;
		for(j = 1; j <= cnt; j++){
			t = i * p[j];
			if(t > maxn) break;
			x[t] = 1;
			if(i % p[j] == 0) break;
			mu[t] = -mu[i];
		}
	}
	for(i = 1; i <= maxn; i++) s[i] = s[i - 1] + mu[i];
	for(i = 1; i <= maxn; i++){
		for(l = 1; l <= i; l = r + 1){
			r = i / (i / l);
			g[i] += (r - l + 1) * (i / l);
		}
	} 
	scanf("%d", &T);
	while(T--){
		scanf("%d%d", &n, &m);
		if(n > m) swap(n, m);
		ans = 0;
		for(l = 1; l <= n; l = r + 1){
			r = min(n / (n / l), m / (m / l));
			ans += z * (s[r] - s[l - 1]) * g[n / l] * g[m / l];
		}
		printf("%lld\n", ans);
	}
	return 0;
}