AtCoder Beginner Contest 162 E - Sum of gcd of Tuples (Hard)

题目描述

一个长度为 \(n\) 的序列 \(\{A_1, A_2, ..., A_n\}\), 其中每一位 \(A_i\) 取值为 \([1, k]\).
共有 \(k^n\) 种序列， 求所有序列的 \(gcd(A_1, A_2, ..., A_n)\) 的和.
\(( 2 \le n \le 10^5, 1 \le k \le 10^5)\)

Solution

\(ans = \sum\limits_{A_1 = 1}^k\sum\limits_{A_2 = 1}^k...\sum\limits_{A_n = 1}^kgcd(A_1, A_2, ..., A_n)\)
\(= \sum\limits_{t = 1}^kt\sum\limits_{A_1 = 1}^k\sum\limits_{A_2 = 1}^k...\sum\limits_{A_n = 1}^k[gcd(A_1, A_2, ..., A_n) = t]\)
\(= \sum\limits_{t = 1}^kt\sum\limits_{A_1 = 1}^{\lfloor\dfrac{k}{t}\rfloor}\sum\limits_{A_2 = 1}^{\lfloor\dfrac{k}{t}\rfloor}...\sum\limits_{A_n = 1}^{\lfloor\dfrac{k}{t}\rfloor}[gcd(A_1, A_2, ..., A_n) = 1]\)
\(= \sum\limits_{t = 1}^kt\sum\limits_{A_1 = 1}^{\lfloor\dfrac{k}{t}\rfloor}\sum\limits_{A_2 = 1}^{\lfloor\dfrac{k}{t}\rfloor}...\sum\limits_{A_n = 1}^{\lfloor\dfrac{k}{t}\rfloor}\sum\limits_{d|gcd(A_1, A_2, ..., A_n)}\mu(d)\)
\(= \sum\limits_{t=1}^kt\sum\limits_{d=1}^{\lfloor\dfrac{k}{t}\rfloor} \mu(d) {\lfloor\dfrac{\lfloor\dfrac{k}{t}\rfloor}{d}\rfloor}^n\)
\(= \sum\limits_{t=1}^kt\sum\limits_{d=1}^{\lfloor\dfrac{k}{t}\rfloor} \mu(d) {\lfloor\dfrac{k}{td}\rfloor}^n\)
\({\lfloor\dfrac{k}{td}\rfloor}\) 和 \({\lfloor\dfrac{k}{t}\rfloor}\) 都可以分块， 复杂度 \(O(n)\)

Code

int n, k;
int primes[N], cnt, mu[N], sum[N];
bool st[N];

void e_prime()
{
	mu[1] = 1;
	for(int i = 2; i < N; ++ i)
	{
		if(!st[i]) 
		{
			primes[cnt ++] = i;
			mu[i] = -1;
		}
		for(int j = 0; i * primes[j] < N; ++ j)
		{
			st[primes[j] * i] = 1;
			if(i % primes[j] == 0) break;
			mu[i * primes[j]] = -mu[i]; 
		}
	}
	for(int i = 1; i < N; ++ i) sum[i] = sum[i - 1] + mu[i];
}

int g(int a, int x)
{
	return a / (a / x);
}

int q_pow(int a, int b)
{
	int res = 1;
	while(b)
	{
		if(b & 1) res = (LL)res * a % P;
		a = (LL)a * a % P;
		b >>= 1;
	}
	return res;
}

int main()
{
	e_prime();
	IOS; cin >> n >> k;
	LL res = 0;
	for(int p = 1, q; p <= k; p = q + 1)
	{
		q = min(k, g(k, p));
		LL A = (LL)(p + q) * (q - p + 1) / 2 % P;
		for(int l = 1, r; l <= k / p; l = r + 1)
		{
			r = min(k / p, g(k / p, l));
			LL B = sum[r] - sum[l - 1];
			LL C = q_pow((k / p) / l, n);
			res = (res + A * B % P * C % P) % P; 
		}
	}
	res = (res % P + P) % P;
	cout << res << endl;
	return 0;
}