【算法讲24：Min25筛】洛谷P5325

内容比(xiang)较(dang)复杂，所以借鉴了很多博文/评论的整合，自己目前只感觉理解了最终 $60\%$ 的内容，毕竟才学了一天半，慢慢加油吧。

引入

• 【6053. 简单的函数】
  给你一个算数函数 $f$，满足：
  $$\begin{array}{rlr}& 1.f(1)=1& \\ & 2.f(p^c)=p\oplus c& \\ & 3.f(ab)=f(a)f(b)& (a\perp b)\end{array}$$
  你需要求出
  $$\sum _{i=1}^{n}f(i)(\mathrm{m}\mathrm{o}\mathrm{d}1e9+7)$$
  数据范围：$1\le n\le 10^{10}$

Min25 筛

• $Min25$ 筛是以亚线性的方式计算出积性函数前缀和。
  特殊条件：$\forall p\in Prime$，$f(p)$ 可以用多项式表示，且 $f(p^c)$ 可以快速得到。
• 我们先考虑一个积性函数 $F(x)$ ，求解：
  $$\sum _{i=1}^n[i\in Prime]F(i)$$
• 设 $P$ 是质数集合， $P_i$ 表示第 $i$ 个质数。
  设：
  $$g(n,j)=\sum _{i=1}^n{i}^k[i\in Pormin(p)>P_j，p|i，p\in P]$$
  翻译成人话：$i$ 是质数，或者 $i$ 的最小质因子大于 $P_j$
  如果我们求的 $f(p)$ 是一个多项式，那么拆成单项式，算出一次项、二次项等的和，最后考虑系数。
• 然后考虑 $g(n,j)$ 函数的转移。
  1. 如果$P_j^2>n$，此时最小质因子是 $P_j$ 的最小合数就是 $P_j^2$
     也就是如果 $P_j^2>n$，此时不会产生新的贡献，有$g(n,j)=g(n,j-1)$
  2. 如果 $P_j^2\le n$，此时 $g(n,j)=g(n,j-1)-X$
     我们减去的贡献就是那些最小质因子是 $P_j$ 的东西
     容斥一下，我们先算上所有含有 $P_j$ 这部分的贡献 $g(n/P_j,j-1)$
     再减掉其他质数以及最小质因数小于 $P_j$ 的那部分，也就是 $g(P_j-1,j-1)$
  3. 最后综合得出转移：
     $$g(n-j)=\{\begin{array}{llc}& g(n,j-1)& P_j^2>n\\ & g(n,j-1)-P_j^k[g(n/P_j,j-1)-g(P_j-1,j-1)]& P_j^2\le n\end{array}$$
     注意这里为了美观， $n/P_j$ 是表示整除的。
  4. 考虑一下 $g(P_j,j)={\sum }_{i=1}^n{P}_j^k$
  5. 我们要求的是 $g(n,|P|)$，其中 $|P|$ 是 $Prime$ 集合的大小。
     根据上面的转移，发现需要的质数只有不大于 $\sqrt{n}$ 的，所以只需要筛出这些质数。
  6. 现在来思考一下 $g$ 函数所代表的含义，可以理解为在模拟埃氏筛的过程
     $g(n,j)$ 表示 $[1,n]$ 排成一列，但是已经筛过一些质数了，把前 $j$ 个质数的倍数全都划掉了，剩下的求个 $F(x)$ 的和就是 $g$ 函数。
     所以转移的过程可以理解为已经晒完了前 $j-1$ 个质数，现在考虑删除第 $j$ 个质数的过程。
  7. 这样枚举量还是很大。我们有：$\lfloor \frac{\lfloor \frac{n}{a}\rfloor }{b}\rfloor =\lfloor \frac{n}{ab}\rfloor$
     所以所有能用到的数都形如 $\lfloor \frac{n}{a}\rfloor$
     用 $mp/ump$ 时间复杂度比较大，我们用手动离散化，对应后面代码里的 $match1[x]为x离散化的下标、match2[x]为\lfloor \frac{n}{x}\rfloor 的下标$
• 那么怎么得到我们要的东西呢？
  1. 设：
     $$S(n,j)=\sum _{i=1}^{n}F(i)[min(p)\ge P_j,p\in P，p|i]$$
     那么，$S(n,j)$ 分为两个部分计算，一部分是所有的质数的和，一部分是所有合数的和，1的值单独计算。
  2. 所有质数的值可以用 $g(n,|P|)$ 表示出来。
  3. 对于合数，枚举一下每个合数的最小质因子以及最小质因子的次幂，进行转移：
     $$S(n,j)=g(n,|P|)-\sum _{i=1}^{j}f(P_i)+\sum _{k\ge j}\sum _{P_k^{e+1}<n}(F(P_k^e)S(\frac{n}{P_k^e},k+1)+F(P_k^{e+1}))$$
     因为 $F(x)$ 是一个积性函数，所以我们把它的最小质因数拆出来，考虑剩下的部分再乘起来。
     因为最小质因数已经被除完，所以剩下部分中不能再含有最小质因数。
     同时所有的 $F(p_k^i)$ 也被筛掉了，需要额外的补进来。
• 代码实现部分
  1. 从 $g(n,j)$ 的递推关系中，我们可以看出：
     首先可以压掉第二维数组
     然后不需要求出所有的 $g_x$，只需要求 $g_{\lfloor \frac{n}{x}\rfloor }$
  2. 但是由于 $N$ 很大，光把下标变成 $\lfloor \frac{n}{x}\rfloor$ 还不够，还要离散化。
     具体的，我们把 $\lfloor \frac{n}{x}\rfloor >\sqrt{n}$ 和 $\lfloor \frac{n}{x}\rfloor \le \sqrt{n}$ 两种分成两类编号来离散化，把空间复杂度压下去。
  3. 由于在 $g$ 的转移方程中我们是把各个质因子拆开来想的，$g_1/g_2$ 代表一次项和二次项。
  4. 有些题目需要特殊考虑2这个唯一的偶质数。
  5. 时间复杂度为 $O(\frac{n^{3/4}}{\log n})$

例子

• P5325 【模板】Min_25筛
  定义积性函数 $f(x)$，且：
  $$f(p^k)=p^k(p^k-1)$$
  求：
  $$\sum _{i=1}^{n}f(i)$$
• 代码：

#include <bits/stdc++.h>
using namespace std;
#define ll long long
const int MAX = 1e6+50;
const int MOD = 1e9+7;
const ll INV6 = 166666668;
ll n,sqrtn;
ll f(ll x){
    x %= MOD;
    return x * (x-1) % MOD;
}
ll pri[MAX],sump1[MAX],sump2[MAX],pcnt;
bool vis[MAX];
void shai(int up){
    for(int i = 2;i <= up;++i){
        if(!vis[i]){
            pri[++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 * pri[j] <= up;++j){
            vis[i * pri[j]] = true;
            if(i % pri[j] == 0)break;
        }
    }
}
ll g1[MAX],g2[MAX],match1[MAX],match2[MAX],w[MAX],tot;
void pre(){
    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] < sqrtn)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;
    }
}
void work(){
    for(int i = 1;i <= pcnt && pri[i] * pri[i] <= n;++i){
        for(int j = 1;j <= tot && pri[i] * pri[i] <= w[j];++j){
            ll tmp2 = w[j] / pri[i];
            ll tmp  = (tmp2 >= sqrtn) ? (match2[n/tmp2]) : (match1[tmp2]);
            g1[j] = (g1[j] - pri[i] * (g1[tmp]-sump1[i-1]+MOD) % MOD  % MOD + MOD) % MOD;
			g2[j] = (g2[j] - pri[i]*pri[i]%MOD * (g2[tmp]-sump2[i-1]+MOD) % MOD + MOD) % MOD;
        }
    }
}
ll S(ll N, ll y){//最后求答案的
	if(pri[y] > N) return 0;
	ll k = (N >= sqrtn) ? (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 * pri[i]*pri[i] <= N; ++i){
		for(ll pk = pri[i]; pri[i]*pk <= N; pk *= pri[i]){
			cpans = (cpans + (f(pk) * S(N/pk, i+1) % MOD + f(pri[i]*pk)) % MOD) %MOD;
		}
	}
	return (prans+cpans) % MOD;
}
int main(){
    scanf("%lld",&n);
    sqrtn = sqrt(n);
    shai(MAX-5);
    pre();
    work();
    cout << (S(n,1) + 1) % MOD;
    return 0;
}

参考

• 小蒟蒻yyb的博客
  Min25筛小结
  洛谷的的题解
  zimindaada、AzusaCat 等用户的题解