Min_25筛学习笔记

发表于 | 分类于

一直想学洲阁筛吧……（这就是我学Min_25筛的理由？）

然后今天终于能yy出来洲阁筛的复杂度证明了……（然后不想写（逃

然后发现了一种更优越的求积性函数前缀和的方法，叫做Min_25筛的……（虽然算是烂大街了？）

Min_25筛

假设现在有一个积性函数\(f\)，我们知道他在质数的幂\(p^c\)上的表现\(f(p^c)\)，然后要去求它的前缀和。

Min_25筛有个重要的要求就是……\(f(p)\)（\(p\)为质数）事一个低次多项式或者一个方便求前缀和的东西。下面就谈谈为什么……

我们先考虑对于所有状态\(\lfloor\frac{n}{x}\rfloor\)，求出范围内所有质数的答案。

看起来蛮棘手的……但我们考虑模仿洲阁筛，定义一个状态\(g(i, x)\)，表示\(x\)范围内事质数或者和前\(i\)个质数都互质的数的答案之和。如果说\(p_i\)达到了\(\lfloor\sqrt{n}\rfloor\)，那么显然\(g(i, x)\)就是\(x\)范围内所有质数的答案了（\(p_i\)表示第\(i\)个质数，下同）。

那么考虑转移（接下来出现的\(f(p_i)\)是指质数的答案多项式中的某一项，因为答案多项式要分开考虑）……对所有\(x\geq p_i^2\)，使用\(g(i, x) = g(i - 1, x) - f(p_i)(g(i - 1,\lfloor\frac{x}{p_i}\rfloor) - g(i - 1, p_i - 1))\)转移（之所以后面又补了一个东西是因为要把多减去的含小于\(p_i\)的质数的合数补回来）；至于\(x < p_i^2\)，我们会发现答案以后就不会变了，因为答案全部都是一堆质数的答案的和（其中不大于\(p_i\)的一定不是合数，否则必定有小于等于\(p_i\)的质因子；倘使大于\(p_i\)的话，没有小于等于\(p_i\)的质因子也只能有且仅有一个大于\(p_i\)的质因子，当然就是质数了）。

这一步操作和洲阁筛几乎如出一辙，复杂度证明也完全一致（复杂度为\(O(\frac{n^{\frac{3}{4}}}{\ln n})\)）。

但是光考虑质数的答案没完惹……再定义状态\(f(i, x)\)表示\(x\)范围内和小于\(p_i\)的质数都互质的数的答案，最终答案显然就是\(f(1) + f(1, n)\)。那么考虑一个很暴力的策略：

先把质数的答案都算进来（之前预处理了），然后考虑合数的答案。考虑枚举所有不小于\(p_i\)的质数\(p\)（假设是第\(j\)个质数），如果说\(p^2 > x\)了那么就没法往下转移了（因为这样用上\(p\)就没法构造合数了），break出来就行了；反之则枚举\(p\)在数中所占的正指数\(e\)，对于所有\(p^{e + 1}\leq x\)，对答案做\(f(j + 1, \lfloor\frac{x}{p^e}\rfloor)f(p^e) + f(p^{e + 1})\)的贡献（其实就是枚举是否选够了\(e\)个\(p\)，然后考虑只用\(p\)的若干次方的情况）。

这个方法看起来极其暴力（甚至也没记忆化），但是复杂度很玄学（也很优秀）……可以参考朱老大的集训队论文。

例题：SPOJ DIVCNT3

这个题显然有\(f(p^c) = 3c + 1\)，然后根据上面说的搞就行了……

代码：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <algorithm>
#include <functional>
#include <utility>
const int maxn = 400005;
using ll = long long;
int prm[maxn]; ll A[maxn], B[maxn];
ll sqr(ll x) {
  ll ret = 1;
  while((ret + 1LL) * (ret + 1LL) <= x) ret ++;
  return ret;
}

ll n, S;
int cnt;
void process() {
  S = sqr(n); cnt = 0;
#ifdef LOCAL
  printf("S : %lld\n", S);
#endif
  for(int i = 1; i <= S; i ++) {
    A[i] = i;
  }
  for(int i = 1; i <= S; i ++) {
    B[i] = n / (ll(i));
  }
  for(int i = 2; i <= S; i ++) {
    if(A[i] == A[i - 1]) continue;
    ll v = A[i - 1], lim = (ll(i)) * (ll(i));
    for(int j = 1; j <= (S / i); j ++) {
      B[j] -= B[j * i] - v;
    }
    for(int j = S / i + 1; j <= S; j ++) {
      ll th = n / (ll(j)); if(th < lim) break;
      B[j] -= A[th / (ll(i))] - v;
    }
    for(int j = S; (ll)j >= lim; j --) {
      A[j] -= A[j / (ll(i))] - v;
    }
    prm[++ cnt] = i;
  }
  prm[++ cnt] = S + 1LL;
}
inline ll query(ll x) {
  if(x <= S) {
    return A[x];
  } else {
    return B[n / x];
  }
}
ll calc(ll m, int x) {
  if(m < (ll)prm[x]) return 0;
  ll ret = 4LL * (query(m) - query(prm[x] - 1));
  for(int i = x; i <= cnt; i ++) {
    ll mul = prm[i];
    if(mul * (ll(prm[i])) > m) break;
    for(int j = 1; ; j ++) {
      if(mul * (ll(prm[i])) > m) break;
      ret += calc(m / mul, i + 1) * (ll(3 * j + 1)) + (ll(3 * j + 4));
      mul *= (ll(prm[i]));
    }
  }
  return ret;
}

int main() {
  int T; scanf("%d", &T);
  while(T --) {
    scanf("%lld", &n);
    process();
    printf("%lld\n", calc(n, 1) + 1LL);
  }
  return 0;
}