Min_25筛学习笔记

Min_25筛可以用来求解一类 \(f(n)\) 的前缀和，其中 \(f\) 是积性函数， \(f(P)\) 是关于质数 \(P\) 的低阶多项式，且 \(\forall k\in N^+,f(P^k)\) 可以通过某些方式快速求解。​

理论部分

首先我们求解 \(\sum_{i=1}^n [i \in \mathbb{P}]f(i)\) 。

首先由于 \(f(P)\) 是多项式，所以我们可以把 \(f(P)\) 的每一项拆开来。

我们令 \(g_{n,j} = \sum_{i=1}^n [i \in \mathbb{P} \space or \space MinFactor(i) > P_j ]i^k\) ，其中 \(MinFactor(i)\) 为 \(i\) 的最小质因子（后简写成 \(MF\)），\(P\) 为质数集，\(P_j\) 为第 \(j\) 个质数。由定义可以得到 \(g_{n,|P|}\) 即为我们所求, \(g_{n,0} = \sum_{i=2}^{n}i^k\)。

那么我们该怎么转移呢？

我们分类考虑一下。如果 \(P_j^2 > n\) ，那么就不存在 \(MF(i) > P_j\) 这条限制了，毕竟不存在一个合数有大于 \(\sqrt n\) 的质因子。

如果 \(P_j^2 \leq n\) ，则由于 \(g_{n,j}\) 所囊括的被计算 \(i^k\) 的 \(i\) 必定是 \(g_{n,j-1}\) 所囊括的的一个子集，所以 \(g_{n,j} = g_{n,j-1}-\Delta\)。

具体考虑 \(\Delta\) 所包含的东西。容易发现 \(g_{n,j}\) 所没有包含的，而 \(g_{n,j-1}\) 却包含的就是最小质因子恰好等于 \(P_j\) 的。

\[∴\Delta = \sum_{i=1}^n [MF(i)=p_j]i^k \\ =\sum_{i=1}^{\lfloor\frac{n}{P_j}\rfloor}[MF(i)\geq p_j] i^k \cdot P_j^k \\=P_j^k\cdot\sum_{i=1}^{\lfloor\frac{n}{P_j}\rfloor}[MF(i)=p_j]i^k\\=P_j^k\cdot(g_{\lfloor\frac{n}{P_j}\rfloor,j-1}-\sum_{i=1}^{j-1}[i\in \mathbb{P}]i^k) \]

把 \(f\) 的那些各个项合在一起，得到转移方程为：

\[g_{n,j}=\begin{cases} g_{n,j-1}&,P_j>\sqrt n\\ g_{n,j-1}-f(P_j)\times(g_{\lfloor \frac{n}{P_j}\rfloor,j-1}-\sum_{i=1}^{j-1} f(P_i))&,P_j\le \sqrt n \end{cases} \]

其中最后的质数前缀和可以在筛质数的时候顺带处理，复杂度不计。

在整完了质数部分的答案后我们就可以合起来，考虑整体的答案了。

我们令\(S_{n,j}=\sum_{i=1}^n[MF(i)>P_j]f(i)\)，则由定义有\(ans = S_{n,0}+f(1)\)，则

\[S_{n,j} = (g_{n,|P|}-\sum_{i=1}^j P_i) + (\sum_{k=j+1}^{P_k\leq \sqrt n} \sum_{e=1}^{P_k^e \leq n}S_{n,k+1} + f(p_k^{e+1})) \\ = (g_{n,|P|}-\sum_{i=1}^j P_i) + (\sum_{k=j+1}^{P_k\leq \sqrt n} \sum_{e=1}^{P_k^e \leq n}f(p_k^e)\cdot S_{\lfloor\frac{n}{p_k^e}\rfloor,k+1} + f(p_k^{e+1})) \]

（第一个式子能倒腾到第二个是因为 \(f(p_k^e)\) 与 \(S_{\lfloor\frac{n}{p_k^e}\rfloor,k+1}\) 互质）

第二个式子第一个括号内的是前文中的质数部分的和，第二个括号内的则是和数部分的和。

第二个括号内的第一个 \(\Sigma\) 相当于在枚举最小质因子，第二个则是在枚举最小质因子的指数。这个 \(S_{\lfloor\frac{n}{p_k^e}\rfloor,k+1}\) 是个递归，由于我们不可能从后面的东西推前面，所以我们枚举他的最小质因子和指数，然后只需要考虑除完之后剩下部分的答案就好了，即 \(\lfloor\frac{n}{p_k^e}\rfloor\) 。因为最小质因数已经被除完，所以剩下部分中不能再含有最小质因数，所以我们递归下去的部分要考虑的最小质因子门槛要提高，即 \(k+1\) 。

同时，为了补上所有被筛掉的 \(p_k^e\) 类的合数，我们要把他们的值加回去，即式子最后的 \(f(p_k^{e+1}))\) 。

然而，我们看到，光是推这个 \(g\) 的时空复杂度，貌似是 \(O(n|P|)\) 的，这很明显无法被接受 比暴力还慢的屑算法有什么存在的必要吗 ，所以说我们在实际的代码中肯定需要一些优化。

代码实现部分

首先我们看到 \(g_{n,j}\) 在 \(j\leq \sqrt n\) 的递推那里，即 \(g_{n,j}=g_{n,j-1}-f(P_j)\times(g_{\lfloor \frac{n}{P_j}\rfloor,j-1}-\sum_{i=1}^{j-1} f(P_i))\) 。我们可以看出，

1. 从 \(j\) 的角度看， \(g_{n,j}\) 都是从 \(g_{x,j-1}\) 转移来的；

2. 从 \(n\) 的角度看， \(g_{n,j}\) 都是从 \(g_{\lfloor \frac{n}{P_j}\rfloor, y}\) 转移来的；

   这意味着什么呢？

从第一点看，我们可以压掉第二维数组。

从第二点看，我们不需要求出所有的 \(g_{x}\) ，只需要求 \(g_{\lfloor \frac{n}{x}\rfloor}\) 。由整除分块的知识得到，这个东西只有可能有 \(\sqrt n\) 种取值。于是时间复杂度降到了 \(O(\sqrt n \times |P|) < O(n)\) （因为预处理的质数只用到 \(\sqrt n\) 就行）。

但是由于 \(n\) 很大，所以我们光把下标变成 \(\lfloor \frac{n}{x}\rfloor\) 不够，还要再离散化。具体的，我们把 \(\lfloor \frac{n}{x}\rfloor > \sqrt n\) 和 \(\lfloor \frac{n}{x}\rfloor \leq \sqrt n\) 两种分成两类编号来“离散化”，这样子就把空间复杂度压下去了。

由于在 \(g\) 的转移方程中我们是把各个质因子拆开来想的，所以在代码里面也得这么做。具体到下面的代码就是\(g1\)和\(g2\)，分别代表一次项和二次项。

尽管洛谷的模板题不需要特殊考虑2这个唯一的偶质数，但是在有些题比如LOJ6053 简单的函数就要考虑到。

时间复杂度为 \(O(n^{1-\epsilon})\) ，一说 \(O(\frac{n^{\frac{3}{4}}}{\log n})\) ，反正我也不会证QAQ

代码（ 洛谷模板（P5325） ）

#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
namespace ztd{
    using namespace std;
    typedef long long ll;
    template<typename T> inline T read(T& t) {//fast read
        t=0;short f=1;char ch=getchar();double d = 0.1;
        while (ch<'0'||ch>'9') {if (ch=='-') f=-f;ch=getchar();}
        while (ch>='0'&&ch<='9') t=t*10+ch-'0',ch=getchar();
        if(ch=='.'){ch=getchar();while(ch<='9'&&ch>='0') t+=d*(ch^48),d*=0.1,ch=getchar();}
        t*=f;
        return t;
    }
}
using namespace ztd;
const int maxn = 1e5+7;
const int mod = 1e9+7;
const ll INV6 = 166666668;//预处理6的逆元
ll n, sqrn;
inline ll f(ll x){//题面的f函数，其定义域只有prime^k
x %= mod;
    return x*(x-1)%mod;
}
ll prime[maxn], sump1[maxn], sump2[maxn], pcnt;
//sump1和sump2分别是质数的前缀和以及二次前缀和(2^2+3^2+5^2……)
bool vis[maxn];
inline void preprocess_prime(int uplimit){//线性筛筛质数
    for(int i = 2; i <= uplimit; ++i){
        if(!vis[i]){
            prime[++pcnt] = i;
            sump1[pcnt] = (sump1[pcnt-1]+i) % mod;
            sump2[pcnt] = (sump2[pcnt-1]+1ll*i*i%mod) % mod;
        } 
        for(int j = 1; j <= pcnt && i*prime[j] <= uplimit; ++j){
            vis[i*prime[j]] = 1;
            if(i % prime[j] == 0) break;
        }
    }
}

ll g1[maxn<<1], g2[maxn<<1], match1[maxn<<1], match2[maxn<<1], w[maxn<<1], tot;
//g1,g2即上文所说， match1和match2和w就是上面所说的离散化了。
inline void preprocess_g(){//预处理g和离散化
    ll tmp;
    for(ll l = 1, r = 0; l <= n; l = r+1){
        r = n/(n/l); ++tot;
        w[tot] = n/l;
        if(w[tot] < sqrn) match1[n/l] = tot;
        else match2[n/(n/l)] = tot;
        tmp = w[tot] % mod;
        g1[tot] = tmp * (tmp+1) / 2 % mod-1;
        g2[tot] = tmp * (tmp+1) % mod * (tmp*2+1) % mod * INV6 % mod - 1;
    }
}
inline void process_g(){//g的递推
    for(int i = 1; i <= pcnt && prime[i]*prime[i] <= n; ++i){
        for(int j = 1; j <= tot && 1ll*prime[i]*prime[i] <= w[j]; ++j){
            ll pmt = w[j]/prime[i];
            ll tmp = (pmt >= sqrn) ? (match2[n/pmt]) : (match1[pmt]);
            g1[j] = (g1[j] - prime[i] * (g1[tmp]-sump1[i-1]+mod) % mod  % mod + mod) % mod;
            g2[j] = (g2[j] - prime[i]*prime[i]%mod * (g2[tmp]-sump2[i-1]+mod) % mod + mod) % mod;
        }
    }
}
ll S(ll N, ll y){//最后求答案的
    if(prime[y] > N) return 0;
    ll k = (N >= sqrn) ? (match2[n/N]) : (match1[N]);
    ll prans = (g2[k]-g1[k]+mod - (sump2[y-1]-sump1[y-1])+mod) % mod;//指质数方面的结果
    ll cpans = 0;//合数方面的结果
    for(int i = y; i <= pcnt && 1ll * prime[i]*prime[i] <= N; ++i){
        for(ll pk = prime[i]; prime[i]*pk <= N; pk *= prime[i]){
            cpans = (cpans + (f(pk) * S(N/pk, i+1) % mod + f(prime[i]*pk)) % mod) %mod;
        }
    }
    return (prans+cpans) % mod;
}
signed main(){
    read(n);  sqrn = sqrt(n);
    preprocess_prime(maxn-5);
    preprocess_g();
    process_g();
    cout << (S(n,1) + 1) % mod << '\n';//这道题定义了f(1)=1
    return 0;
}