（赛前练手 #9） BZOJ2005 [Noi2010]能量采集（容斥原理）

2005: [Noi2010]能量采集
Time Limit: 10 Sec Memory Limit: 512 MB
Submit: 5041 Solved: 3078
[Submit][Status][Discuss]
Description
栋栋有一块长方形的地，他在地上种了一种能量植物，这种植物可以采集太阳光的能量。在这些植物采集能量后，
栋栋再使用一个能量汇集机器把这些植物采集到的能量汇集到一起。 栋栋的植物种得非常整齐，一共有n列，每列
有m棵，植物的横竖间距都一样，因此对于每一棵植物，栋栋可以用一个坐标(x, y)来表示，其中x的范围是1至n，
表示是在第x列，y的范围是1至m，表示是在第x列的第y棵。 由于能量汇集机器较大，不便移动，栋栋将它放在了
一个角上，坐标正好是(0, 0)。 能量汇集机器在汇集的过程中有一定的能量损失。如果一棵植物与能量汇集机器
连接而成的线段上有k棵植物，则能量的损失为2k + 1。例如，当能量汇集机器收集坐标为(2, 4)的植物时，由于
连接线段上存在一棵植物(1, 2)，会产生3的能量损失。注意，如果一棵植物与能量汇集机器连接的线段上没有植
物，则能量损失为1。现在要计算总的能量损失。 下面给出了一个能量采集的例子，其中n = 5，m = 4，一共有20
棵植物，在每棵植物上标明了能量汇集机器收集它的能量时产生的能量损失。 在这个例子中，总共产生了36的能
量损失。
Input
仅包含一行，为两个整数n和m。

Output
仅包含一个整数，表示总共产生的能量损失。

Sample Input
【样例输入1】

5 4

【样例输入2】

3 4

Sample Output
【样例输出1】

36

【样例输出2】

20

对于100%的数据：1 ≤ n, m ≤ 100,000。

HINT
Source

第一道NOI题目AC祭…（看来是比较简单的题目了）
首先我们需看清问题求解的本，本题我们可以很轻松地发现，每个点对答案的贡献恰好是2 * gcd(i , j) - 1，y，此时问题就转换了求${\sum }_{i=1}^n{\sum }_{j=1}^{m}gcd(i,j)$
做法有二：
①莫比乌斯反演 + 整除分块。做的时候考虑过但还是不会推。。所以请大家还是自己搜一下反演的做法吧
②容斥原理。曾经做过一道n = m的题，于是把那道题的方法套用过来，其实对于gcd(x , y) = k来说(x <= n && y <= m)，一定也有gcd(x / k , y / k) = 1，那么我们就可以先枚举x的因数，也就是k，然后运用质因数分解和容斥原理来筛有多少对数满足它们的gcd为k
AC Code:

#include<iostream>
#include<cstdio>
#include<cstring>
#include<cmath>
#define rg register
#define il inline
#define maxn 200005
#define ll long long
using namespace std;
bool vis[maxn];
int primes[maxn] , f[maxn] , phi[maxn] , qq[maxn] , q[maxn] , cnt , tot , tot2;
ll n , m;
il ll read(){
	rg ll x = 0 , w = 1 ;
	rg char ch = getchar();
	while(ch < '0' || ch > '9'){
		if (ch == '-') w = -1;
		ch = getchar();	
	}
	while (ch >= '0' && ch <= '9'){
		x = (x << 3) + (x << 1) + ch - '0';
		ch = getchar();	
	}
	return x * w;
}
void prime(int n){
	for (rg int i = 2 ; i <= n ; ++i){
		if (!vis[i]) primes[++cnt] = i , phi[i] = i - 1;
		for (rg int j = 1 ; j <= cnt && i * primes[j] <= n ; ++j){
			vis[i * primes[j]] = 1;
			if (!(i % primes[j])){
				phi[i * primes[j]] = phi[i] * primes[j];
				break;
			}
			phi[i * primes[j]] = phi[i] * (primes[j] - 1);
		}
	}
}
void fenjie(int x){
	tot = 0;
	rg int tmp = x , i = 1;
	while (vis[tmp]){
		if (!(tmp % primes[i])){
			q[++tot] = primes[i];
			tmp /= primes[i];	
		}
		while (!(tmp % primes[i]))
			tmp /= primes[i];
		++i;
	}
	if (tmp != 1){
		q[++tot] = tmp;	
	}
}
void fenjie2(int x){
	tot2 = 0;
	rg int sqx = sqrt(x);
	for (rg int i = 1 ; i <= sqx ; ++i){ 
		if (!(x % i)){
			qq[++tot2] = i , qq[++tot2] = x / i;
			if (i == x / i) --tot2;	
		}
	}
}
int rongchi(int x){
	f[0] = -1;
	rg int ans = 0;
	cnt = 0;
	rg int tmp = 0;
	for (rg int i = 1 ; i <= tot ; ++i){
		tmp = cnt;
		for (rg int j = 0 ; j <= tmp ; ++j){
			f[++cnt] = f[j] * q[i] * (-1);
		}
	}
	for (rg int i = 1 ; i <= cnt ; ++i)
		ans += x / f[i];
	return x - ans;
}
int gcd(int a , int b){
	if (!b) return a;
	return gcd(b , a % b);	
}
int min(int a , int b){
  	return (a < b) ? a : b;
}
int max(int a , int b){
  	return (a > b) ? a : b;
}
int main(){
	n = read() , m = read();
	prime(max(n , m));
	if (n > m) n ^= m , m ^= n , n ^= m;
	ll ans = m , tmp;
	for (rg int i = 2 ; i <= n ; ++i){
		tmp = 0;
		if (!vis[i]) {
			ans += m / i * i + m - m / i;
			continue;
		}
		fenjie2(i);
		for (rg int j = 1 ; j <= tot2 ; ++j){
			fenjie(i / qq[j]);
			tmp = rongchi(m / qq[j]);
			ans += qq[j] * tmp;
		}
	}	
	printf("%lld" , ans * 2 - n * m);
	return 0;
}