min_25筛学习笔记+牛客多校02E

本来没有学习这种较难的算法的想法的，因为比赛也做不到这种难度的题， 但是最近打牛客多校02，有一题要求 $[1,n]$ 中素数的个数，我以为是像莫反一样容斥，但是后面感觉不行。赛后知道是用 min_25 筛来求，赛时过了一车，因此我也不得不学习这个算法了。

我打算拿洛谷里面的模板来举例。

就是给你积性函数
$$f(x)=\{\begin{array}{ll}k& x=1\\ g(x)={\sum }_{i=0}^k{a}_i{x}^i& x=p^c,p\in prime,c\ge 1\\ f(x_1)f(x_2)& (x_1,x_2)=1\end{array}$$

然后让你求
$$\sum _{i=1}^{n}f(i)$$

其中 $n\le 10^{13}$

在牛客多校的那个题中，也有 $n\le 10^9$

一般来讲，想解决这种问题都要从 $i$ 最小的质因子入手，因为一个数要么是质数，要么最小的质因子 $\le n^{0.5}$，而 $n^{0.5}$ 一般都是$10^6$~$10^7$ 数量级，可以直接用线性筛得到，也能全部枚举一次。因此总的思路大概就是先求出 ${\sum }_{i=1}^{n}g(i)$ 然后通过枚举最小的质因子来修正有超过一个质因子的数的贡献，这样能够很好的解决大质数无法直接得到的问题。

如果没有学过要怎么求 ${\sum }_{i=1}^n{i}^k$，可以看我的博客

基于这个大的思路，我自己感觉 min_25 筛的思路是非常自然，很容易就能记下来。

对于式子
$$\sum _{i=1}^{n}f(i)$$

先将 $1$ 去掉，因为其不包含任何素数，计算的时候容易出现奇怪的错误，虽然按照积性函数的定义， $f(1)=1$。

于是变成
$$f(1)+\sum _{i=2}^{n}f(i)$$

自己感觉 min_25 筛很关键的一点是其能够快速求出
$$g(b)=\sum _{p\in prime,p\le b}{p}^k$$其中 $b=\lfloor \frac{n}{a}\rfloor$，这也意味着 $b$ 的个数只有 $O(n^{0.5})$

虽然只能求出若干次方的和，但只要把每一项都求一次就可以得到原函数了。

然后剩下的和数就能通过枚举小最小的质因子来考虑

因此 min_25 筛有两步，分别是求出质数出的函数和，求出所有的函数和

第一步：求出质数处的函数和

这一步的思路就是像我刚才说的，先求出所有的函数和，然后通过枚举最小的质因子把和数的贡献剪掉即可。

我们从小到大枚举素数，并且同时维护所有的 $g(b)$，假设当前的素数为 $p$，假设其为第 $i$ 个素数。那么显然有 $p\le n^{0.5}$。 我们记第 $i$ 个素数为 $p_i$

对于 $g(b)$，我们应该减去 $1$ ~ $b$ 中最小质因子为 $p$ 的和数，如果我们直接固定 $p$，也就是只需要考虑 $1$ ~ $\lfloor \frac{b}{p}\rfloor$ 中最小质因子 $\ge p$ 的数，然后你就会发现这个数刚好就是 $g(\lfloor \frac{b}{p}\rfloor )-{\sum }_{j=1}^{i-1}{p}_j^k$，（因为 $g$ 函数中包含 $<p$ 的素数）后者只需要线性筛的时候预处理即可。

然后容斥完以后我们就可以得到所有的 $g(b)$ 了。

第二步：求出原函数

这一步的思路是也是拆分成质数贡献与和数贡献，质数贡献可以直接通过前面的 $g$ 函数得到，和数贡献可以通过枚举最小质因子不断递归求解。

我们另 $S(n,x)$ 表示 $1$ ~ $n$ 中最小质因子 $>p_x$ 的数的贡献之和，我们认为质数本身就是其最小质因子。

然后有
$$S(n,x)=g(n)-\sum _{i=1}^{x}f(p_i)+\sum _{p_i^c\le n,i>x}f(p_i^c)(S(\lfloor \frac{n}{p_i^c}\rfloor ,i+1)+[c>1])$$

因为 $c>1$ 的时候 $p_i^c$ 本身也能够产生贡献。

以下是我洛谷模板的代码

#include <bits/stdc++.h>
#define ll long long
using namespace std;const int mod = 1e9 + 7;int power(int a, int b) {int ret = 1;while (b) {if (b & 1) ret = 1ll * ret * a % mod;a = 1ll * a * a % mod; b >>= 1;}return ret;
}const int inv2 = power(2, mod - 2), inv3 = power(3, mod - 2);
const int N = 1e6 + 10;int cnt, p[N]; bool v[N];
int sp1[N], sp2[N];
ll n, Sqr, w[N];
int tot, g1[N], g2[N], ind1[N], ind2[N];void pre(int n) {v[1] = 1;for (int i = 2; i <= n; i++) {if (!v[i]) {p[++cnt] = i;sp1[cnt] = (sp1[cnt - 1] + i) % mod;sp2[cnt] = (sp2[cnt - 1] + 1ll * i * i) % mod;}for (int j = 1; j <= cnt && i * p[j] <= n; j++) {v[i * p[j]] = 1;if (i % p[j] == 0) break;}}
}int S(ll x, int y) {if (p[y] >= x) return 0;ll k = x <= Sqr ? ind1[x] : ind2[n / x];int ans = (1ll * g2[k] - g1[k] + mod - (sp2[y] - sp1[y]) + mod) % mod;for (int i = y + 1; i <= cnt && 1ll * p[i] * p[i] <= x; i++) {ll pe = p[i];for (int e = 1; pe <= x; e++, pe *= p[i]) {ll xx = pe % mod;ans = (ans + xx * (xx - 1) % mod * (S(x / pe, i) + (e != 1)) % mod) % mod;}}return ans % mod;
}int main() {// freopen("in.in", "r", stdin);ios::sync_with_stdio(false);cin.tie(nullptr);cin >> n;Sqr = sqrt((long double)n);pre(Sqr);for (ll i = 1, j; i <= n; i = j + 1) {j = n / (n / i);w[++tot] = n / i;int now = w[tot] % mod;g1[tot] = 1ll * now * (now + 1) / 2 % mod - 1;g2[tot] = 1ll * now * (now + 1) / 2 % mod * (2 * now + 1) % mod * inv3 % mod - 1;if (n / i <= Sqr) ind1[n / i] = tot;else ind2[n / (n / i)] = tot;}for (int i = 1; i <= cnt; i++) {for (int j = 1; j <= tot && 1ll * p[i] * p[i] <= w[j]; j++) {ll k = w[j] / p[i] <= Sqr ? ind1[w[j] / p[i]] : ind2[n / (w[j] / p[i])];g1[j] = (g1[j] - 1ll * p[i] * (g1[k] - sp1[i - 1] + mod) % mod + mod) % mod;g2[j] = (g2[j] - 1ll * p[i] * p[i] % mod * (g2[k] - sp2[i - 1] + mod) % mod + mod) % mod;}}cout << (S(n, 0) + 1) % mod << "\n";return 0;
}

时间复杂度

第一步的时间复杂度是 $O(\frac{n^{0.75}}{\ln n})$ 是没有争议的，有一些人底下写的是 $\log$ 也是对的，因为只有常数倍的差别

对于第一重循环，因为素数的密度的是 $\frac{1}{\ln n}$，所以可以提出一个 $\frac{1}{\ln n}$

然后变为这个式子 $O({\sum }_{i=1}^{\sqrt{n}}\sqrt{i}+\sqrt{\frac{n}{i}})=O({\sum }_{i=1}^{\sqrt{n}}\sqrt{\frac{n}{i}})=O(\sqrt{n}{\sum }_{i=1}^{\sqrt{n}}\sqrt{\frac{1}{i}})=O(\sqrt{n}{\int }_1^n\frac{1}{\sqrt{x}}dx)=O(\sqrt{n}\sqrt{n}{|}_1^{\sqrt{n}})=O(n^{0.75})$

第二步我在网上没有找到具体的证明，自己也感觉比较难证明，但是听说是一般情况是 $O(n^{1-\delta })$，但是如果 $n\le 10^{13}$，就满足 $O(\frac{n^{0.75}}{\ln n})$，基本适合OI的范围。

学完 min_25 筛，就会发现前面那个求素数个数的问题实在是太简单了，连 $S$ 函数都不用算，直接把 $g$ 函数计算以下就行了，然后把多项式设置为 $f(x)=1$ 即可。