题解 P4450 【双亲数】

· · 题解

upd：256 天后的更新，修复了错误和公式格式，充实了内容

如果 \LaTeX 渲染挂了请移步洛谷博客或个人博客

under 题解 P4450

蒟蒻第一道莫反题

本题解适合新手学习基础，dalao 请跳过，不要嘲讽 qwq

前置知识：

单位函数（\epsilon）：

\epsilon(i) = \begin{cases} 1 &, i = 1 \\ 0 &, i \ne 1 \end{cases}

即当且仅当 i = 1 时取值取值为 1，否则为 0。很简单对不对？

Dirichlet 卷积(\ast)：

Dirichlet 卷积是数论函数之间的运算。即

h = f \ast g

其中 f, g, h 均为数论函数。则

h(n) = \sum_{d | n}{f(d) g \left(\frac nd \right)}

也可以写作

h(n) = \sum_{i · j = n}{f(i) g(j)}

我们发现：\epsilon 函数是 Dirichlet 卷积的单位元，即 f = f \ast \epsilon。

Mobius 函数(\mu)：

\mu(n) = \begin{cases} 1 &, n = 1 \\ (-1)^m &, n = p_1p_2 \cdots p_m \\ 0 &, \mathrm{otherwise} \end{cases}

其中 p_1,p_2, \cdots, p_m 是不同的质数。\mu(n) 恰在 n 无平方因子时取值非零。显然 \mu 是积性函数。

\begin{aligned} \sum_{d | n}{\mu(d)} &= \epsilon(n) \\ \epsilon &= 1 \ast \mu \end{aligned}

Mobius 变换

设 f 是数论函数，若数论函数 g 满足

\begin{aligned} g &= f * \mathbf1 \\ \textrm{That is} \quad g(n) &= \sum_{d | n}{f(d)} \end{aligned}

则称 g 是 f 的 Mobius 变换。

推一下式子。

g = f \ast 1 \Rightarrow f = f \ast \epsilon = f \ast 1 \ast \mu = g \ast \mu

即

g = f \ast 1 \Leftrightarrow f = g \ast \mu

至此就可以来做这道题了。

回到了原题。题目要求的是

\sum_{i = 1}^A{\sum_{j = 1}^B{[\gcd(i, j) = d]}}

对式子进行套路性的变形

\sum_{i = 1}^{\lfloor A / d \rfloor}{\sum_{j = 1}^{\lfloor B / d \rfloor}{[\gcd(i, j) = 1]}}

出现了 [ \gcd(i, j) = 1]，也即 \epsilon( \gcd(i, j))。把 \epsilon = 1 \ast \mu 代入，得到了

\sum_{i = 1}^{\lfloor A / d \rfloor}{\sum_{j = 1}^{\lfloor B / d \rfloor}{\sum_{g | \gcd(i, j)}{\mu(g)}}}

（似乎越来越丑了）此时可以把 d | \gcd(i, j) 拆开：

\begin{aligned} & \sum_{i = 1}^{\lfloor A / d \rfloor}{\sum_{j = 1}^{\lfloor B / d \rfloor}{\sum_{g | i, g | j}{\mu(g)}}} \\ =& \sum_{g = 1}^{\min(A / d, B / d)}{\mu(g) \left\lfloor \frac{A}{dg} \right\rfloor \left\lfloor \frac{B}{dg} \right\rfloor} \end{aligned}

至此，计算答案可以 O(n) 解决。

算法还可以继续优化。发现式子里还有两个形如 \lfloor n / d \rfloor 的因子，因此可以整除分块优化成 O(\sqrt{n}) 的。配合 O(n) 的线性筛 \mu 函数，总复杂度 O(n)，可以轻松通过本题。

code 1

std::bitset<1000003> v;
int p[78501], mu[1000003];

void sieve(int N) {
  int m = 0;

  v.set(1), mu[1] = 1;
  for (int i = 2; i <= N; ++i) {
    if (!v[i]) p[++m] = i, mu[i] = -1;
    for (int j = 1; j <= m && i * p[j] <= N; ++j) {
      v.set(i * p[j]);
      if (i % p[j] == 0) break;
      mu[i * p[j]] = -mu[i];
    }
    mu[i] += mu[i - 1];
  }
}

int main() {
  int A, B, d;
  long long ans = 0;

  std::cin >> A >> B >> d;
  A /= d, B /= d;
  if (A > B) std::swap(A, B);

  sieve(A);
  for (int l = 1, r, a, b; l <= A; l = r + 1) {
    a = A / l, b = B / l;
    r = std::min(A / a, B / b);
    ans += 1ll * a * b * (mu[r] - mu[l - 1]);
  }
  std::cout << ans << '\n';
  return 0;
}

code 2

如果你像我一样闲得难受，可以学来杜教筛写写这题，复杂度好像还更高了。

不会杜教筛但想学的朋友可以来看我的这篇题解。

这个多次询问的题目让我对杜教筛及空间优化有了更深入的理解（

std::bitset<10003> v, vm;
int MX, base, p[1231], mu[10003], mu_[10003];

inline int to(int x) {
  return MX / x;
}

void sieve(int N); // 同 code 1

int getMu(int N) {
  if (N <= base) return mu[N];
  if (vm[to(N)]) return mu_[to(N)];
  int ans = 0;
  for (int l = 2, r, L; l <= N; l = r + 1) {
    L = N / l, r = N / L;
    ans += (r - l + 1) * getMu(L);
  }
  vm.set(to(N));
  return mu_[to(N)] = 1 - ans;
}

int main() {
  int A, B, d;
  long long ans = 0;

  std::cin >> A >> B >> d;
  A /= d, B /= d;
  if (A > B) std::swap(A, B);

  MX = A, base = pow(A, 2. / 3);
  sieve(base);
  for (int l = 1, r, a, b, gR, gL; l <= A; l = r + 1) {
    a = A / l, b = B / l;
    r = std::min(A / a, B / b);
    vm.reset(), gR = getMu(r);
    vm.reset(), gL = getMu(l - 1);
    ans += 1ll * a * b * (gR - gL);
  }
  std::cout << ans << '\n';
  return 0;
}