Codeforces Round 911 (Div. 2) D

D. Small GCD

题意

定义\(f(a,b,c)\)为\(a,b,c\)中较小两个数的\(gcd\),给定数组\(a_{1...n}\),求\(\sum\limits_{i=1}^{n}\sum\limits_{j=i+1}^{n}\sum\limits_{k=j+1}^{n}f(a_i,a_j,a_k)\)

题解

显然可以先排序,不会影响结果,排完序后\(a_k\)就是最大的,不会对\(gcd(a_i,a_j)\)产生影响.

所以当我们去枚举中间的\(a_j\)时,那么对于\(a_j\)来说,产生\((n-j)\times \sum\limits_{i = 1}^{n-1}gcd(a_i,a_j)\),其中\((n-j)\)是因为后面的\((n-j)\)个\(a_k\)都不会对\(gcd(a_i,a_j)\)产生影响,所以答案最终就是\(\sum\limits_{i=1}^{n-2}\sum\limits_{j=i+1}^{n-1}gcd(a_i,a_j) \times (n-j)\),而这样做的复杂度是\(O(n^2logm)\),其中\(m\)是\(\max\limits_{1\leq i\leq n}\{a_i\}\),而这样是会超时的.

到了这里,就该推出今天刚了解的欧拉反演了!

欧拉反演:即\(n\)的所有因子的欧拉函数和为n.

\(\sum\limits_{d|n}\varphi(d)=n\).

将\(n\)换成其他:

\(gcd(i,j) = \sum\limits_{d|gcd(i,j)}\varphi(d) = \sum\limits_{d|i}\sum\limits_{d|j}\varphi(d)\)

则:

\(\sum\limits_{i=1}^{n}gcd(i,n) = \sum\limits_{i=1}^{n}\sum\limits_{d|i}\sum\limits_{d|n}\varphi(d)=\sum\limits_{d|n}\sum\limits_{i=1}^{n}\sum\limits_{d|i}\varphi(d)=\sum\limits_{d|n}\frac{n}{d}\varphi(d)\)

即:

\(\sum\limits_{i=1}^{n}gcd(i,n)=\sum\limits_{d|n}\frac{n}{d}\varphi(d)\)

因为\(1\sim n\)里面含有因子为\(d\)一共有\(\frac{n}{d}\)个,所以这里就直接替换了,不过在这题里面从\(a_1 \sim a_{j-1}\)里面并不知道含有\(d\)作为因子数的有多少,所以我们需要维护\(a_1 \sim a_{j-1}\)中每个数的所有因子的个数\(cnt_d\),那么要计算\(1 \sim j-1\)的\(gcd\)则可以替换成:

\(\sum\limits_{i=1}^{j-1}gcd(a_i,a_j) = \sum\limits_{d|a_j}cnt_d\varphi(d)\)

对于每一个数,先预处理出所含的约数,然后维护\(cnt_d\),就可以直接利用公式求和了.

ACcode:

#include <bits/stdc++.h>
#define debug(a) cout<<#a<<"="<<a<<'\n';

using namespace std;
using i64 = long long;

typedef pair<i64, i64> PII;

vector<int> euler_range(int n) {
	vector<int> phi(n + 1), prime;
	vector<bool> is_prime(n + 1, true);
	is_prime[1] = 0, phi[1] = 1;
	for (int i = 2; i <= n; i++) {
		if (is_prime[i]) prime.push_back(i), phi[i] = i - 1;
		for (int j = 0; j < (int)prime.size() && i * prime[j] <= n; j++) {
			is_prime[i * prime[j]] = 0;
			if (i % prime[j]) phi[i * prime[j]] = phi[i] * (prime[j] - 1);
			else {
				phi[i * prime[j]] = phi[i] * prime[j];
				break;
			}
		}
	}
	return phi;
};

constexpr int N = 1E5;

auto phi = euler_range(N);
vector<int> fac[N + 1];

void solve() {

	int n;
	cin >> n;
	vector<i64> a(n), b;
	for (auto &i : a) cin >> i;
	ranges::sort(a);
	i64 ans = 0;
	vector<int> cnt(N + 1);
	for (int j = 0; j < n; j ++) {
		for (auto d : fac[a[j]])
			ans += 1ll * phi[d] * (n - j - 1) * cnt[d];
		for (auto d : fac[a[j]])
			cnt[d] ++;
	}

	cout << ans << '\n';

}
int main() {
	ios::sync_with_stdio(false);
	cin.tie(nullptr);

	for (int i = 1; i <= N; i ++)
		for (int j = 1; j <= N / i; j ++)
			fac[i * j].push_back(i);

	int T;
	cin >> T;
	while (T--) {
		solve();
	}

	return 0;
}