首页 > 其他分享 >「学习笔记」Min_25 筛

「学习笔记」Min_25 筛

时间:2023-01-18 21:45:48浏览次数:49  
标签:25 le Min int sum 笔记 MAXN minp mathrm

学长来讲 Min_25 筛了。

然后我咕了快一周了,终于写了!!

前情提要:没有复杂度证明。

Min_25 筛干啥的?

问 Min_25 啊

Min_25 筛是用来在 \(O(\frac{n^{\frac{3}{4}}}{\log n})\) 的复杂度内解决一些积性函数前缀和的问题。

先摆出这个积性函数的要求:

\(p\) 是质数,那么要求有:

  1. \(f(p)\) 是一个关于 \(p\) 的多项式
  2. \(f(p^k)\) 便于计算。

符号规定

  • \(p_i\) 代表第 \(i\) 个质数
  • \(\mathrm{minp}(i)\) 代表 \(i\) 的最小质因子

Part I

首先,我们定义一个函数 \(S(x,j)\):

\[S(x,j)=\sum_{i=2}^xf(i)[\mathrm{minp}(i)\ge p_j] \]

既:\([2,x]\) 中最小质因子大于等于第 \(j\) 个质数的所有数的函数值的和。

我们先不考虑 \(1\),那么发现答案就是 \(S(x,1) + 1\)。

然后把这东西展开:

\[\begin{aligned} S(x,j)&=\sum_{i=2}^xf(i)[\mathrm{minp}(i)\ge p_j]\\ &= \underbrace{\sum_{\substack{i\ge j\\p_i\le n}} f(p_i)}_{\textrm{所有质数处的答案}} + \underbrace{\sum_{\substack{k\ge j\\p_k^2\le n}}\sum_{\substack{e\ge 1\\p_k^{e+1}\le n}}}_{\textrm{枚举最小的质因子与指数,即 $p_k^e$}}\underbrace{f(p_k^e)S\left(\left\lfloor\frac{x}{p_k^e}\right\rfloor, k + 1\right)}_{\textrm{所有形如 $p_k^e\times \cdots$} 的数} + \underbrace{f(p_k^{e+1})}_{\textrm{所有形如 $p_k^{e+1}$} 的数} \end{aligned}\]

其实就是把所有数分成了这三种类型:

  1. \(p_k\):这个计算方法后面会讲。
  2. \(p_k^e\times \cdots\):由于是积性函数,直接乘即可。剩下的数最小质因子就必须大于 \(p_k\) 了,所以是 \(S\left(\left\lfloor\dfrac{x}{p_k^e}\right\rfloor, k + 1\right)\)。
  3. \(p_k^e (e\ge 2)\):直接计算即可。式子中用 \(e+1\) 的原因是防止重复计算 \(e=1\) 的情况。

然后大概讲一下中间两个 \(\sum\) 的限制条件:

  1. \(p_k^2 \le n\):如果 \(p_k^2 > n\),那么首先 \(p_k^{e+1} > n\),这一项就不存在了;并且 \(p_k \times \cdots > p_k^2 > n\),所以第一项也不存在,所以只需要计算 \(p_k^2 \le n\) 的情况就可以了。
  2. \(p_k^{e+1}\le n\):当 \(p_k^{e+1}> n\) 时,形如 \(p_k^e \times \cdots\) 的数也一定大于 \(n\),所以当 \(p_k^{e+1}\le n\) 时,两个式子都会卡到上界。

于是这一大堆式子就是 \(O(\frac{n^{\frac{3}{4}}}{\log n})\) 的。

我不会证

Part II

那么如何计算 \(\sum_{i\ge j,p_i\le n} f(p_i)\)?

首先我们可以转换一下,改成前缀和相减的形式。\(j\) 比较小,\(\sum_{\substack{i<j}} f(p_i)\) 是可以预处理出来的,所以我们只需要求出 \(\sum_{i\ge 1,p_i\le n} f(p_i)\) 即可。

我们再设一个函数 \(g(x,j)\):

\[g_k(x,j)=\sum_{i=2}^xi^k[i \in \mathrm{prime}\ \mathrm{or}\ \mathrm{minp}(i)>p_j] \]

那么发现:当 \(j\rightarrow +\infty\) 时,就不存在 \(i\) 满足 \(\mathrm{minp}(i)>p_j\) 了,也就是:

\[g_k(x,+\infty)=\sum_{i=2}^xi^k[i \in \mathrm{prime}] \]

由于 \(f(p)\) 是一个多项式,那么我们只需要把若干个 \(g_k(x,+\infty)\) 加起来,就是我们要求的答案。

我们可以考虑递归地求这个东西。

当 \(p_j^k > n\) 时,发现此时就不存在满足 \(\mathrm{minp}(i)>p_j\) 且不是质数的数了,所以: $$g_k(x,j)=g_k(x,j-1)$$

当 \(p_j^k \le n\) 时,我们考虑 \(\mathrm{minp}(i)>p_j\) 在 \(\mathrm{minp}(i)>p_{j-1}\) 的基础上多筛掉了哪些数。首先肯定 \(\mathrm{minp}(i)=p_j\),那么因为 \(i^k\) 是完全积性函数,所以我们直接乘上 \(p_j^k\),然后得到这样的式子:

\[g_k(x,j)=g_k(x,j-1) - p_j^k\left(g_k\left(\left\lfloor\frac{x}{p_j}\right\rfloor,j-1\right)-\sum_{i=1}^{j-1}p_i^k\right) \]

因为 \(p_j\) 的指数不一定只有 \(1\),所以这里就体现出完全积性函数的好处了:可以只提出来一个 \(p_j\),后面可以随意乘,这样得到的数一定有 \(p_j\) 而且指数不一定。但是我们需要减去所有质数处的值,因为前后的质数都需要计算。

于是这个式子可以直接进行递推求,然后我们再带回原来的式子就好了。

Part III

代码实现:

#include <bits/stdc++.h>
using namespace std;
const int MAXN = 1000005, P = 1000000007;
long long n;
int sqr;
int pri[MAXN], pcnt, m;
bool vis[MAXN];
long long w[MAXN];
int g1[MAXN], g2[MAXN];
int p1[MAXN], p2[MAXN];
int mp[2][MAXN];
int qpow(int a, int b) {
    int ans = 1;
    while (b) {
        if (b & 1) ans = 1ll * ans * a % P;
        a = 1ll * a * a % P;
        b >>= 1;
    }
    return ans;
}
void seive(int n) {
    for (int i = 2; i <= n; i++) {
        if (!vis[i]) {
            pri[++pcnt] = i;
        }
        for (int j = 1; j <= pcnt && i * pri[j] <= n; j++) {
            vis[i * pri[j]] = 1;
            if (i % pri[j] == 0) break;
        }
    }
    for (int i = 1; i <= pcnt; i++) p1[i] = (p1[i - 1] + pri[i]) % P;
    for (int i = 1; i <= pcnt; i++) p2[i] = (p2[i - 1] + 1ll * pri[i] * pri[i] % P) % P;
}
int getId(long long x) { // 离散化
    return x <= sqr ? mp[0][x] : mp[1][n / x];
}
const int INV2 = qpow(2, P - 2), INV6 = qpow(6, P - 2);
void getG(long long n) {
    for (long long l = 1, r; l <= n; l = r + 1) { // 离散化
        r = n / (n / l);
        w[++m] = n / l;
        if (w[m] <= sqr) mp[0][w[m]] = m;
        else mp[1][n / w[m]] = m;
        g1[m] = (w[m] % P * ((w[m] + 1) % P) % P * INV2 % P + P - 1) % P;
        g2[m] = (w[m] % P * ((w[m] + 1) % P) % P * ((2 * w[m] % P + 1) % P) % P * INV6 % P + P - 1) % P;
    }
    for (int j = 1; j <= pcnt; j++) {
        for (int i = 1; i <= m && 1ll * pri[j] * pri[j] <= w[i]; i++) {
            int id = getId(w[i] / pri[j]);
            g1[i] = (g1[i] - 1ll * pri[j] * (g1[id] - p1[j - 1] + P) % P + P) % P;
            g2[i] = (g2[i] - 1ll * pri[j] * pri[j] % P * (g2[id] - p2[j - 1] + P) % P + P) % P;
        }
    }
}
int S(long long x, int j) {
    if (x <= 1 || pri[j] > x) return 0;
    int id = getId(x), ans = ((1ll * (g2[id] - g1[id]) - (p2[j - 1] - p1[j - 1])) % P + P) % P;
    for (int k = j; k <= pcnt && 1ll * pri[k] * pri[k] <= x; k++) {
        long long p1 = pri[k], p2 = p1 * p1;
        for (int e = 1; p2 <= x; e++, p1 = p2, p2 *= pri[k]) {
            ans = (ans + p1 % P * (p1 % P - 1) % P * S(x / p1, k + 1) % P + p2 % P * (p2 % P - 1) % P) % P;
        }
    }
    return ans;
}
int main() {
    scanf("%lld", &n);
    sqr = sqrt(n);
    seive(sqr);
    getG(n);
    printf("%d\n", (S(n, 1) + 1) % P);
    return 0;
}

标签:25,le,Min,int,sum,笔记,MAXN,minp,mathrm
From: https://www.cnblogs.com/apjifengc/p/17060635.html

相关文章