[LibreOJ 2185][SDOI2015]约数个数和

danihao123 posted @ 2018年7月10日 21:23 in 题解 with tags SDOI 莫比乌斯反演 狄利克雷卷积 loj , 509 阅读
转载请注明出处:http://danihao123.is-programmer.com/

首先我们有这样一个式子:

\[d(ij) = \sum_{a | i}\sum_{b | j} [\gcd(\frac{i}{a}, b) = 1]\]

怎样理解这一结论呢?我们知道\(ij\)的约数\(d\)一定可以被分解为\(ab\)使得\(a | i, b | j\),但这种划分可能会非常的多,我们要保证每个约数\(d\)只有一种划分被统计过。

那么考虑上述划分方法。对于\(i\)和\(j\)中一个有一个没有的质因子,是一定会划分到\(a\)和\(b\)中固定的一个的。那么考虑一个在\(i\)和\(j\)中都出现的质因子\(p\),假设\(p\)在\(d\)中出现的次数超过了\(i\)和\(j\)中任何一个的话,那么\(p\)在\(a\)中一定会尽可能多的出现(否则\(\frac{i}{a}\)中会出现\(p\),而这样无论如何一定有\(p | b\),因此会不满足\(\gcd(\frac{i}{a}, b) = 1\))。如果没超过呢?那么一定会全部划分到\(a\)中(其他情况下有\(p | b\)并且会出现\(p | \frac{i}{a}\))。

下面考虑应用这个结论(为了方便,使用等价的\(d(ij) = \sum_{a | i}\sum_{b | j} [\gcd(a, b) = 1]\),假设\(n\le m\)):

\[
\begin{aligned}
\quad&\sum_{i = 1}^n\sum_{j = 1}^m \sum_{a | i}\sum_{b | j}\sum_{d | i, d | j}\mu(d)\\
=&\sum_{d = 1}^n\mu(d)\sum_{a = 1}^{\lfloor\frac{n}{d}\rfloor}\lfloor\frac{n}{ad}\rfloor\sum_{b = 1}^{\lfloor\frac{m}{d}\rfloor}\lfloor\frac{m}{bd}\rfloor\\
=&\sum_{d = 1}^n\mu(d)\sum_{a = 1}^{\lfloor\frac{n}{d}\rfloor}d(a)\sum_{b = 1}^{\lfloor\frac{m}{d}\rfloor}d(b)
\end{aligned}
\]

后面两个和式可以随便\(O(n)\)欲处理一下(然后我当初写了个\(O(n\sqrt{n})\)的神必预处理……),然后这个题做完力,,,

代码:

#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cctype>
#include <algorithm>
#include <queue>
#include <utility>
typedef long long int_t;
const int maxn = 50005;
int_t miu[maxn], miu_s[maxn], f[maxn];
void process() {
  static bool vis[maxn];
  static int prm[maxn]; int cnt = 0;
  const int N = 50000;
  miu[1] = 1LL;
  for(int i = 2; i <= N; i ++) {
    if(!vis[i]) {
      miu[i] = -1LL; prm[cnt ++] = i;
    }
    for(int j = 0; j < cnt && prm[j] * i <= N; j ++) {
      int v = prm[j] * i;
      vis[v] = true;
      if(i % prm[j] == 0) {
        miu[v] = 0; break;
      } else {
        miu[v] = miu[i] * -1;
      }
    }
  }
  for(int i = 1; i <= N; i ++) {
    miu_s[i] = miu_s[i - 1] + miu[i];
  }
  for(int i = 1; i <= N; i ++) {
    for(int j = 1; j <= i;) {
      int nx = i / (i / j);
      f[i] += (int_t(nx - j + 1)) * (int_t(i / j));
      j = nx + 1;
    }
  }
}

int_t solve(int n, int m) {
  int_t ans = 0;
  for(int i = 1; i <= std::min(n, m);) {
    int nx = std::min(m / (m / i), n / (n / i));
    ans += (miu_s[nx] - miu_s[i - 1]) * f[m / i] * f[n / i];
    i = nx + 1;
  }
  return ans;
}
#ifdef LOCAL
#define LO "%I64d"
#else
#define LO "%lld"
#endif
int main() {
  process(); int T;
  scanf("%d", &T);
  while(T --) {
    int n, m; scanf("%d%d", &n, &m);
    printf(LO, solve(n, m)); puts("");
  }
  return 0;
}

登录 *


loading captcha image...
(输入验证码)
or Ctrl+Enter