HDU 5528 Count a * b 欧拉函数

题意：

定义函数\(f(n)\)为\(i \cdot j \not\equiv 0 \; (mod \; n)\)的数对\((i,j)\)的个数\((0 \leq i,j \leq n)\)
\(g(n)=\sum_{d|n}f(d)\)，求\(g(n) \; mod \; 2^{64}\)，其中\(1 \leq n \leq 10^9\)

仰慕叉老师手推公式

分析：

计算\(f(n)\)

\(f(n)\)不容易计算，但可以计算它的反面，也就是计算满足\(ij \equiv 0 \; (mod \; n)\)的数对\((i,j)\)的个数。
为了方便计算，对\(i,j\)的范围稍作改动，\(1 \leq i,j \leq n\)，这不会影响最终计算结果。
\(i \cdot j\)是\(n\)的整数倍，也就是\(n\)被分到了这里面。
不妨设\(i=ax,j=by\)，其中\(ab=n\)。我们可以枚举\(a\)来统计答案，但是这样有重复计算。
比如\((3 \cdot 2) \cdot (2 \cdot 1) = (6 \cdot 1) \cdot (2 \cdot 1) = 6 \cdot 2\)，这样数对\((6, 2)\)在\(a=3\)和\(a=6\)时被计算了两次。
所以计数时我们规定，数对\((i,j)\)只会在\(a\)最大的时候被计算，所以有\(gcd(b,x)=1\)成立。
答案为\(\sum\limits_{a|n} \sum\limits_{1 \leq x \leq \frac{n}{a}}[gcd(x,b)=1] \cdot \frac{n}{b}=\sum\limits_{a|n}a \phi(\frac{n}{a})\)

计算\(g(n)\)

根据题中给的关系，得到\(g(n)=\sum\limits_{m|n}(m^2-f(m))=\sum\limits_{m|n}(m^2-\sum\limits_{d|m}d \phi(\frac{m}{d}))\)
继续化简后面那部分，\(\sum\limits_{m|n}\sum\limits_{d|m}d \phi(\frac{m}{d})\)交换求和顺序，得到\(\sum\limits_{d|n}d\sum\limits_{\frac{m}{d}|\frac{n}{d}}\phi(\frac{m}{d})\)
由于\(\sum\limits_{d|n} \phi(d)=n\)，所以得到\(\sum\limits_{d|n}d \cdot \frac{n}{d}=n\sum\limits_{d|n}1=n \sigma(n)\)。其中\(\sigma(n)\)为\(n\)的约数的个数。

\(\sigma(n)\)可以这样计算：
首先将\(n\)质因数分解为\(\prod\limits_{i=1}^{r}p_{i}^{e_i}\)，这样\(\sigma(n)=\prod\limits_{i=1}^{r}(e_i+1)\)。
因为我们已经将\(n\)分解，所以计算前面那部分的时候不要在\(\sqrt{n}\)内枚举\(n\)的约数，会超时的。要根据分解出来的结果直接枚举。

#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
using namespace std;

const int maxn = 32000 + 10;
int pcnt, prime[maxn];
bool vis[maxn];

void preprocess() {
	for(int i = 2; i < maxn; i++) {
		if(!vis[i]) prime[pcnt++] = i;
		for(int j = 0; j < pcnt && i * prime[j] < maxn; j++) {
			vis[i * prime[j]] = true;
			if(i % prime[j] == 0) break;
		}
	}
}

int tot, p[maxn], e[maxn];

unsigned long long ans;
void dfs(int d, unsigned long long product) {
	if(d == tot) {
		ans += product * product;
		return ;
	}

	dfs(d + 1, product);
	for(int i = 1; i <= e[d]; i++) {
		product *= p[d];
		dfs(d + 1, product);
	}
}

int main()
{
	preprocess();
	int T; scanf("%d", &T);
	while(T--) {
		int n;
		scanf("%d", &n);

		int t = n;
		tot = 0;
		for(int i = 0; i < pcnt && prime[i] * prime[i] <= n; i++) {
			if(t % prime[i] == 0) {
				p[tot] = prime[i];
				e[tot] = 0;
				while(t % prime[i] == 0) {
					t /= prime[i];
					e[tot]++;
				}
				tot++;
			}
		}
		if(t > 1) { p[tot] = t; e[tot++] = 1; }

		unsigned long long tmp = 1;
		for(int i = 0; i < tot; i++) tmp = tmp * (e[i] + 1);
		tmp *= n;

		ans = 0;
		dfs(0, 1);

		ans -= tmp;

		printf("%llu\n", ans);
	}

	return 0;
}