Codeforces Gym103428(2021CCPC威海) I. Distance
Codeforces Gym103428 I. Distance
分析
贪心地思考，从 \(i\) 走向 \(j\) ，我们可以先从 \(i\) 走向 \(\gcd(i,j)\) ，再从 \(\gcd(i,j)\) 走向 \(j\) 。
方法一
从 \(i\) 走向 \(\gcd(i,j)\) 的过程中，中间会经过乘积为 \(\frac{i}{\gcd(i,j)}\) 的若干边（其中每条边的边权为素数），实际上走过的距离就是将 \(\frac{i}{\gcd(i,j)}\) 进行质因数分解，然后将每个质因子（含次数）求和。我们不妨把某个数\(x\)进行如上过程记为函数 \(f(x)\) 。
正式地，\(f(x)\) 的定义如下：
\[ x=p_1^{k_1} p_2^{k_2}\cdots p_t^{k_t}\\
f(x)=k_1p_1+k_2p_2+\cdots+k_tp_t
\]
函数 \(f(x)\) 有一个不错的性质（如下式），利用这个性质可以对目标式子进行化简。
\[ \forall x,y\in \N^* ,f(xy)=f(x)+f(y)
\]
对目标式子的化简如下：
\[ \begin{aligned}
&\sum\limits_{i=1}^n
\sum\limits_{j=1}^n
\left(
f\left(\frac{i}{\gcd(i,j)}\right)
+f\left(\frac{j}{\gcd(i,j)}\right)
\right)\\
=&
\sum\limits_{i=1}^n
\sum\limits_{j=1}^n
\left(
f\left(i\right)
+f\left(j\right)-2f\left(\gcd(i,j)\right)
\right)\\
=&
2n\sum\limits_{i=1}^nf(i)-2\sum\limits_{d=1}^n
f(d)
\sum\limits_{i=1}^{\frac{n}{d}}
\sum\limits_{j=1}^{\frac{n}{d}}
[\gcd(i,j)=1]\\
=&
2n\sum\limits_{i=1}^nf(i)-2\sum\limits_{d=1}^n
f(d)
\left(
2\sum\limits_{i=1}^{\frac{n}{d}}\varphi(i)-1
\right)
\end{aligned}
\]
尽管 \(f(i)\) 不是积性函数，但仍然可以利用洲阁筛预处理 \(f(i)\) 的前缀和的根号个点值，只是需要修改对应的转移方式。
对于上式的后半部分，利用整除分块，配合洲阁筛预处理出的 \(\varphi(i)\) 的前缀和，进行计算。
下面解释洲阁筛中后半部分计算 \(f(i)\) 的转移方式。
记 \(e\) 为质因子次数，\(\mathbf P_i\) 为第 \(i\) 个质数，原本对于积性函数的转移如下：
\[ S(n, j) = S(n, j + 1) + \sum_{\mathbf P_j^{e + 1} \leq n} f(\mathbf P_j^e) \left(S\left( \frac n {\mathbf P_j^e}, j + 1 \right) - S(\mathbf P_j, j + 1)\right) + f(\mathbf P_j^{e + 1})
\]
在此题中，\(f(i)\) 由积性函数的乘法变成了加法，故上面和式中\(f\)与\((S-S)\)相乘的部分需要改成加号，但需要注意的是，不应加一个\(f\)，而应该加 \((S-S)\) 和式中所代表元素的个数个 \(f\)。
所以需要额外计算\((S-S)\) 和式中所代表元素的个数，实际上可以引入一个新函数 \(H\)，其中 \(\mathrm{minp}\) 为最小质因子。
\[ H(n, j) = \sum_{i=1}^n [\mathrm{minp}(i) \geq \mathbf P_j \lor i \in \mathbf P]\cdot 1
\]
显然 \(1\) 是积性函数，可以直接套用上面对于积性函数的转移。
于是\((S-S)\) 和式中所代表元素的个数可以被表示为
\[ H\left( \frac n {\mathbf P_j^e}, j + 1 \right) - H(\mathbf P_j, j + 1)
\]
故转移变成
\[ S(n, j) = S(n, j + 1) + f(\mathbf P_j^e)\left(H\left( \frac n {\mathbf P_j^e}, j + 1 \right) - H(\mathbf P_j, j + 1)\right) +\sum_{\mathbf P_j^{e + 1} \leq n} \left(S\left( \frac n {\mathbf P_j^e}, j + 1 \right) - S(\mathbf P_j, j + 1)\right) + f(\mathbf P_j^{e + 1})
\]
此方式常数有点大，时限4s，将将通过。
方法二
对于 \(\text{dis}(i,j)\) ，我们可以对质因子进行枚举，根据上述贪心的思想，可以发现，对于某个质因子 \(p\) ，其贡献为 \(i, j\) 分别包含质因子 \(p\) 的次数之差的绝对值的 \(p\) 倍。
正式地，记素数集为 \(\mathbb{P}\) ，整数 \(i\) 包含质因子的次数为 \(e_i\) ，上面的想法可以表述为下式：
\[ \text{dis}(i,j) = \sum\limits_{p \in \mathbb{P}}p|e_i - e_j|
\]
于是，容易想到对于目标式子，我们也可以枚举素数，对于一个素数 \(p\) ，其对目标式子的贡献为
\[ 2p\sum_{c=1}^{\lfloor\log_pn\rfloor} \left\lfloor\frac{n}{p^c}\right\rfloor \left(n-\left\lfloor\frac{n}{p^c}\right\rfloor\right)
\]
于是，对于目标式子，可以写成：
\[ 2\sum\limits_{p \in \mathbb{P}}p\sum_{c=1}^{\lfloor\log_pn\rfloor} \left\lfloor\frac{n}{p^c}\right\rfloor \left(n-\left\lfloor\frac{n}{p^c}\right\rfloor\right)
\]
对于小于等于 \(\sqrt{n}\) 的素数，直接按式子暴力计算。
对于大于 \(\sqrt{n}\) 的素数，有 \(c=1\)，直接使用Min_25筛的前半部分计算素数的前缀和的根号个点值，利用整除分块进行计算。
代码
方法一
#include <cmath>
#include <iostream>
using namespace std;
typedef long long Lint;
const int mod = 1e9 + 7;
const int inv2 = 500000004;
inline int add(int a, int b) {
int res = a + b;
return res < mod ? res : res - mod;
}
inline int mul(int a, int b) {
return (Lint)a * b % mod;
}
inline int neg(int a) {
return mod - a;
}
namespace Min_25 {
const int maxm = 1000005;
Lint n;
int m;
int vis[maxm], pri[maxm], sp0[maxm], sp1[maxm], tot;
void euler() {
vis[1] = 1;
for (int i = 2; i <= m; i++) {
if (!vis[i])
pri[++tot] = i;
for (int j = 1; j <= tot && (Lint)i * pri[j] <= m; j++) {
vis[i * pri[j]] = 1;
if (i % pri[j] == 0)
break;
}
}
for (int i = 1; i <= tot; i++) {
sp0[i] = add(sp0[i - 1], 1);
sp1[i] = add(sp1[i - 1], pri[i]);
}
}
Lint w[maxm];
int tot_w;
int g0[maxm], g1[maxm], id1[maxm], id2[maxm];
int f[maxm], phi[maxm], cnt[maxm];
int id(Lint x) {
return x <= m ? id1[x] : id2[n / x];
}
void solve() {
for (Lint l = 1, r, val; l <= n; l = r + 1) {
val = n / l;
r = n / val;
w[++tot_w] = val;
int tmp = w[tot_w] % mod;
g0[tot_w] = add(tmp, neg(1));
g1[tot_w] = add(mul(mul(tmp, add(tmp, 1)), inv2), neg(1));
if (w[tot_w] <= m) {
id1[w[tot_w]] = tot_w;
} else {
id2[n / w[tot_w]] = tot_w;
}
}
for (int i = 1; i <= tot; i++) {
for (int j = 1; j <= tot_w && (Lint)pri[i] * pri[i] <= w[j]; j++) {
Lint num = w[j] / pri[i];
int k = id(num);
g0[j] = add(g0[j], neg(add(g0[k], neg(sp0[i - 1]))));
g1[j] = add(g1[j], mul(pri[i], neg(add(g1[k], neg(sp1[i - 1])))));
}
}
for (int i = 1; i <= tot_w; i++) {
f[i] = g1[i];
cnt[i] = g0[i];
phi[i] = add(g1[i], neg(g0[i]));
}
for (int i = tot; i >= 1; i--) {
for (int j = 1; j <= tot_w && (Lint)pri[i] * pri[i] <= w[j]; j++) {
for (Lint k = pri[i], c = 1, k_phi = pri[i] - 1, k_f = pri[i]; k * pri[i] <= w[j]; k *= pri[i], k_phi = mul(k_phi, pri[i]), k_f = add(k_f, pri[i]), ++c) {
f[j] = add(add(add(f[j], add(f[id(w[j] / k)], neg(f[id(pri[i])]))), add(k_f, pri[i])), mul(k_f, add(cnt[id(w[j] / k)], neg(cnt[id(pri[i])]))));
cnt[j] = add(add(cnt[j], add(cnt[id(w[j] / k)], neg(cnt[id(pri[i])]))), 1);
phi[j] = add(add(phi[j], mul(k_phi, add(phi[id(w[j] / k)], neg(phi[id(pri[i])])))), mul(k_phi, pri[i]));
}
}
}
}
int get_f(Lint n) {
if (n < 2)
return 0;
return f[id(n)];
}
int get_phi(Lint n) {
if (n < 2)
return 1;
return phi[id(n)] + 1;
}
} // namespace Min_25
int main() {
Lint n;
cin >> n;
Min_25::n = n;
Min_25::m = sqrt(n);
Min_25::euler();
Min_25::solve();
int n_mod = n % mod;
Lint ans = mul(n_mod, Min_25::get_f(n));
for (Lint i = 1, r, val; i <= n; i = r + 1) {
val = n / i;
r = n / val;
int sum_f = add(Min_25::get_f(r), neg(Min_25::get_f(i - 1)));
int phi = Min_25::get_phi(val);
ans = add(ans, mul(neg(sum_f), add(phi, add(phi, neg(1)))));
}
cout << add(ans, ans) << '\n';
return 0;
}
方法二
#include <cmath>
#include <iostream>
using namespace std;
typedef long long Lint;
const int maxsqrn = 1000005;
const int mod = 1e9 + 7;
const int inv2 = 500000004;
int vis[maxsqrn], pri[maxsqrn], tot;
int sp[maxsqrn], id1[maxsqrn], id2[maxsqrn], g[maxsqrn];
Lint w[maxsqrn];
int tot_w;
inline int add(int a, int b) {
int res = a + b;
return res < mod ? res : res - mod;
}
inline int mul(int a, int b) {
return (Lint)a * b % mod;
}
inline int neg(int a) {
return mod - a;
}
void euler(int m) {
vis[1] = 1;
for (int i = 2; i <= m; i++) {
if (!vis[i])
pri[++tot] = i;
for (int j = 1; j <= tot && (Lint)i * pri[j] <= m; j++) {
vis[i * pri[j]] = 1;
if (i % pri[j] == 0)
break;
}
}
for (int i = 1; i <= tot; i++) {
sp[i] = add(sp[i - 1], pri[i]);
}
}
void getG(int m, Lint n) {
for (Lint l = 1, r, val; l <= n; l = r + 1) {
val = n / l;
r = n / val;
w[++tot_w] = val;
int tmp = w[tot_w] % mod;
g[tot_w] = add(mul(mul(tmp, add(tmp, 1)), inv2), neg(1));
if (w[tot_w] <= m) {
id1[w[tot_w]] = tot_w;
} else {
id2[n / w[tot_w]] = tot_w;
}
}
for (int i = 1; i <= tot; i++) {
for (int j = 1; j <= tot_w && (Lint)pri[i] * pri[i] <= w[j]; j++) {
Lint num = w[j] / pri[i];
int k = (num <= m) ? id1[num] : id2[n / num];
g[j] = add(g[j], mul(pri[i], neg(add(g[k], neg(sp[i - 1])))));
}
}
}
int main() {
Lint n;
cin >> n;
int m = sqrt(n);
euler(m);
getG(m, n);
int ans = 0;
for (int i = 1; i <= tot; i++) {
int tmp = 0;
for (Lint j = pri[i]; j <= n; j *= pri[i]) {
tmp = add(tmp, mul(n / j % mod, (n - n / j) % mod));
}
ans = add(ans, mul(tmp, pri[i]));
}
for (Lint i = m + 1, pre = sp[tot] % mod, val, vv; i <= n; i = vv + 1) {
val = n / i;
vv = n / val;
int k = (vv <= m) ? id1[vv] : id2[n / vv];
ans = add(ans, mul(add(g[k], neg(pre)), mul(val % mod, (n - val) % mod)));
pre = g[k];
}
cout << mul(ans, 2) << '\n';
return 0;
}
浙公网安备 33010602011771号