说起积性函数求和我最早是在一个在线评测系统里被一道n 10^10的模板题按在地上摩擦的。那时候评论区还有人问“为什么不线性筛到 n 再求前缀和”我差点就信了。等你真把线性筛开到10^7都会卡内存的时候才会理解这类问题的本质求和范围和预处理能力之间存在几个数量级的落差必须用筛法把复杂度压到亚线性。今天要聊的“一种筛法”并不是指筛素数的那套经典思路而是数论里处理大范围积性函数前缀和的两类代表技巧杜教筛和 Min_25 筛。前者适合函数本身能构造出好算的狄利克雷卷积后者适合函数在素数处的取值是多项式形式、但整体又不那么好构造的场景。这两种方法都能把Σ_{i1}^n f(i)的复杂度压到O(n^(2/3))或O(n^(3/4) / log n)的量级让n 10^12也不再是遥不可及的范围。这篇文章面向的是已经会线性筛、懂一点整除分块但对“大范围求和”还没有形成系统思路的读者。我会把原理推导、代码模板、以及我实际踩过的坑一次讲清楚。1. 积性函数求和难在哪1.1 这类问题到底是什么简单说积性函数就是满足f(ab) f(a) * f(b)、且gcd(a, b) 1的数论函数。常见的欧拉函数φ(n)、莫比乌斯函数μ(n)、除数函数σ_k(n)全是积性函数。求和问题就是给定一个很大的n求S(n) Σ_{i1}^{n} f(i)直接枚举 i、对每个数分解质因数再算 f(i)单次分解就是O(√i)整体做下来是O(n√n)哪怕n 10^6都够呛。线性筛可以做到在O(n)内预处理出所有 f(i)但这里有个前提你得真的开一个长度为 n 的数组。竞赛里这类题给的 n 一般是10^9到10^12线性筛连数组都开不下更别提跑完。于是“筛法”在这里的含义就变了我们不再把每个数都筛一遍而是只保留那些“被整除分块整除出来的关键位置”在这些位置上做递推或分治。1.2 线性筛为什么救不了场很多人刚接触这类题时会陷入一个误区既然φ、μ都有线性筛写法那是不是加个前缀和数组就完事了理论上确实如此但代价是空间和时间都是O(n)。当n 10^10时光是long long数组就要80GB这已经超出了绝大多数评测环境的承受范围。更微妙的是即使空间够线性筛的过程也要求从小到大遍历每个数这个遍历本身就有10^10次操作在当前硬件条件下无论如何也是秒级以上的和题目限时完全冲突。所以大范围求和问题几乎是被逼着走亚线性路线核心思想只有一个把所有⌊n/i⌋的取值当作“关键点”在这些点上做递归或迭代把无效的中间值全部跳过去。杜教筛和 Min_25 筛正是沿着这个思路演化出来的两个流派。2. 杜教筛拆卷积换前缀和2.1 从狄利克雷卷积推递推式杜教筛的核心逻辑是给目标函数f找一个辅助积性函数g让h f * g狄利克雷卷积的前缀和特别好算。狄利克雷卷积的定义是h(n) Σ_{d|n} f(d) * g(n/d)然后我们在求和号上做文章Σ_{i1}^{n} h(i) Σ_{d1}^{n} g(d) * S_f(⌊n/d⌋)把这个等式右边第一项d 1那一项单独拎出来就能得到S_f(n) Σ_{i1}^{n} h(i) - Σ_{d2}^{n} g(d) * S_f(⌊n/d⌋)注意右边的Σ_{d2}^{n}可以通过整除分块变成O(√n)个不同的⌊n/d⌋每个递归调用S_f(⌊n/d⌋)再用记忆化保存结果。这样总复杂度就是O(n^(2/3))前提是Σ h(i)能在O(1)或O(√n)内快速计算Σ g(d)在整除分块的每一段前缀和都能快速得到。2.2 三个经典组合直接套实际做题时最常用的组合其实少得可怜我把它们的卷积关系整理成了表目标函数 f辅助函数 gh f * gΣ h(i) 的快速公式μ(n)1常函数ε(n)1φ(n)1常函数id(n) nn(n1)/2n·φ(n)id(n)ε(n)1为什么这么选因为Σ_{i1}^n 1 n和Σ_{i1}^n i n(n1)/2都是闭式而ε的前缀和更是恒等于 1。比如求S_φ(n)时S_φ(n) n(n1)/2 - Σ_{d2}^{n} S_φ(⌊n/d⌋)在这个式子里g(d)恒为 1整除分块每一段的g前缀和就是段长(r - l 1)实现起来非常干净。2.3 实现细节与代码模板写杜教筛时有几个细节直接影响正确性和性能预处理阈值建议取n^(2/3)左右先线性筛出这些范围内的f前缀和超过阈值再用记忆化递归。预处理太少会导致递归层数太多太多则预处理本身消耗过大。记忆化键不要直接用unordered_maplong long, long long存所有值因为递归访问的值只有O(√n)个但每次哈希都有开销。我一般用小范围数组加快表混合处理。整除分块时注意r n / (n / l)这行代码是核心写错会直接死循环。下面是一份能直接跑的杜教筛模板同时求μ和φ的前缀和#include bits/stdc.h using namespace std; using ll long long; const int MAXN 5000000; // 预处理阈值约 n^(2/3) ll phi[MAXN 5], mu[MAXN 5]; ll prefPhi[MAXN 5], prefMu[MAXN 5]; vectorint primes; bool isComp[MAXN 5]; void init() { phi[1] 1; mu[1] 1; for (int i 2; i MAXN; i) { if (!isComp[i]) { primes.push_back(i); phi[i] i - 1; mu[i] -1; } for (int p : primes) { if (1LL * i * p MAXN) break; isComp[i * p] true; if (i % p 0) { phi[i * p] phi[i] * p; mu[i * p] 0; break; } else { phi[i * p] phi[i] * (p - 1); mu[i * p] -mu[i]; } } } for (int i 1; i MAXN; i) { prefPhi[i] prefPhi[i - 1] phi[i]; prefMu[i] prefMu[i - 1] mu[i]; } } unordered_mapll, ll mpPhi, mpMu; ll S_phi(ll n) { if (n MAXN) return prefPhi[n]; if (mpPhi.count(n)) return mpPhi[n]; ll ans n * (n 1) / 2; for (ll l 2, r; l n; l r 1) { r n / (n / l); ans - (r - l 1) * S_phi(n / l); } return mpPhi[n] ans; } ll S_mu(ll n) { if (n MAXN) return prefMu[n]; if (mpMu.count(n)) return mpMu[n]; ll ans 1; // Σ ε(i) 1 for (ll l 2, r; l n; l r 1) { r n / (n / l); ans - (r - l 1) * S_mu(n / l); } return mpMu[n] ans; } int main() { init(); ll n; while (cin n) { cout S_phi(n) S_mu(n) \n; } return 0; }这份模板我用了很久唯一要提醒的是n * (n 1) / 2在n 10^12时会超过long long的舒适区边缘如果在模意义下做必须先把n和n 1分别取模再乘上2的逆元。后面第 5 章会专门展开这个坑。3. Min_25 筛把质数挑出来再拼合数杜教筛很优雅但它的命门在于你得能找到一个让Σ h和Σ g都好算的辅助函数g。如果f(p)在素数处的取值不是一个简单多项式比如f(p) p^2 3p 1这种或者你根本找不到合适的卷积组合杜教筛就卡壳了。这时候就该 Min_25 筛上场。3.1 第一步 G所有 x 处的质数贡献Min_25 筛的思路可以理解成“筛掉合数只留质数”。它先对所有整除块值x ⌊n/i⌋求出 Σ_{p ≤ x} f(p)也就是所有不超过 x 的素数处的函数值之和。做法很巧妙先假装所有数都是质数给每个整除块值初始化一个“全量贡献”然后从小到大枚举素数p把最小质因子恰好是p的合数贡献减掉。每减一轮剩下的就是“还没被更小质数筛掉的数的贡献”最终所有合数的贡献都被去掉剩下的就是真正的质数和。这一步的转移公式长这样G(x, j) G(x, j-1) - f(p_j) * ( G(⌊x/p_j⌋, j-1) - Σ_{ij} f(p_i) )其中j表示已经用前j个质数筛过。由于f(p)是多项式可以拆成若干p^k分别维护每个k维护一组g_k(x) Σ_{i2}^{x} i^k最后按系数组合回去。3.2 第二步 S拼接 p^e 与剩余部分质数部分的贡献拿到之后还差合数部分。Min_25 筛的第二段用递归来枚举“最小质因子及其指数”核心公式是S(x, y) G(x) - sp[y-1] Σ_{k≥y, p_k^2 ≤ x} Σ_{e≥1, p_k^(e1) ≤ x} ( f(p_k^e) * S(⌊x/p_k^e⌋, k1) f(p_k^(e1)) )其中sp[y-1]是前y-1个质数的 f 值之和G(x) - sp[y-1]代表所有不小于第y个质数的质数贡献。合数部分则是枚举最小质因子p_k的指数e把p_k^e提出来剩余部分继续递归约束“最小质因子要大于 p_k”。初看这公式有点吓人但它本质就是一个唯一分解每个合数都能写成“最小质因子的幂次 × 剩余部分”而剩余部分的最小质因子必定更大于是递归天然不重不漏。3.3 实现细节与代码模板写 Min_25 筛最容易出错的地方有两个。一是整除分块值w[t]到数组下标的映射当x ≤ √n时用id1[x]否则用id2[n/x]。二是第一段更新g时对同一个质数p_j必须从大到小遍历整除块防止本轮更新污染后续引用的旧值。这一点和 01 背包倒序遍历是同一个道理。下面是一份求Σ φ(i) mod MOD的 Min_25 筛模板我把f(p) p - 1拆成了g1和g0两个数组分别维护#include bits/stdc.h using namespace std; using ll long long; const int MOD 1000000007; const ll INV2 500000004; ll n, sqrtN; vectorll primes; ll sp1[100005], sp0[100005]; ll w[200005]; int id1[100005], id2[100005]; ll g1[200005], g0[200005]; int tot 0; void init_primes(ll N) { vectorbool comp(N 5, false); for (ll i 2; i N; i) { if (!comp[i]) primes.push_back(i); for (ll p : primes) { if (i * p N) break; comp[i * p] true; if (i % p 0) break; } } for (int i 1; i (int)primes.size(); i) { sp1[i] (sp1[i - 1] primes[i - 1]) % MOD; sp0[i] (sp0[i - 1] 1) % MOD; } } int get_id(ll x) { return x sqrtN ? id1[x] : id2[n / x]; } void prepare_g() { sqrtN sqrt(n); for (ll l 1, r; l n; l r 1) { r n / (n / l); ll x n / l; w[tot] x; if (x sqrtN) id1[x] tot; else id2[n / x] tot; g1[tot] (x % MOD * ((x 1) % MOD) % MOD * INV2 % MOD - 1 MOD) % MOD; g0[tot] (x % MOD - 1 MOD) % MOD; } for (int j 0; j (int)primes.size(); j) { ll p primes[j]; if (p * p n) break; for (int t tot; t 1 w[t] p * p; --t) { int id get_id(w[t] / p); g1[t] (g1[t] - p * ((g1[id] - sp1[j] MOD) % MOD) % MOD MOD) % MOD; g0[t] (g0[t] - (g0[id] - sp0[j] MOD) % MOD MOD) % MOD; } } } ll S(ll x, int y) { if (x 2 || y (int)primes.size()) return 0; if (x primes[y]) return 0; ll res (g1[get_id(x)] - g0[get_id(x)] MOD) % MOD; res (res - ((sp1[y] - sp0[y]) % MOD) MOD) % MOD; for (int k y; k (int)primes.size(); k) { ll p primes[k]; if (p * p x) break; ll pe p; for (int e 1; pe * p x; e, pe * p) { ll f_pe (pe % MOD - (pe / p) % MOD MOD) % MOD; res (res f_pe * S(x / pe, k 1) % MOD) % MOD; res (res (pe % MOD * (p % MOD) % MOD - pe % MOD MOD) % MOD) % MOD; } } return res; } int main() { cin n; init_primes(sqrtN); prepare_g(); ll ans (S(n, 0) 1) % MOD; // 最后补上 f(1) 1 cout ans \n; return 0; }注意这份模板的函数S(x, y)里y用的是 0-based 索引也就是从第y个质数开始考虑。外部调用S(n, 0)表示不考虑任何质数下界的限制最后把1加回去就是完整的Σ φ(i)。4. 怎么选复杂度、边界与实战判断4.1 两种算法的复杂度对比杜教筛的复杂度依赖预处理范围和记忆化的命中率理论上是O(n^(2/3))实际跑n 10^12时预处理到10^8左右比较合适但内存又会成为瓶颈。我的经验是当n不超过10^11且能轻松找到卷积组合时杜教筛的代码更短、调试成本更低。Min_25 筛的时间复杂度约为O(n^(3/4) / log n)从渐进角度不如杜教筛好看但它的常数通常很小而且对函数的要求更宽松。实际跑n 10^10到10^12时Min_25 筛往往比杜教筛更快因为它的第一段预处理把质数部分一次性算清了第二段递归的层数也不深。维度杜教筛Min_25 筛适用范围能找到辅助函数 g使 h 前缀和好算f(p) 是低次多项式f(p^e) 可快速计算时间复杂度O(n^(2/3))依赖预处理阈值O(n^(3/4) / log n)实现难度较低中等典型函数μ、φ、n·φ(n)σ_k、φ、以及各类自定义积性函数记忆化必须否则指数爆炸不一定需要递归天然剪枝4.2 快速判断表与示例数据判断用哪个我一般先做三步如果题目函数是 μ 或 φ且我一眼能写出辅助函数 g直接用杜教筛。如果函数是f(p) p^k的线性组合或题目明确给了f(p^e)的算式直接用 Min_25 筛。如果两个都行选自己代码更熟的那个因为二分的常数和边界处理在赛场上更宝贵。举个例子求Σ_{i1}^{n} σ_2(i)其中σ_2(p^e) 1 p^2 p^4 ... p^(2e)。这个函数在素数处的取值是p^2 1配合常数项可以拆分杜教筛想找辅助函数反而绕用 Min_25 筛就是模板题级别的难度。反过来如果只要求Σ φ(i)杜教筛三行递推式就结束了Min_25 筛反而杀鸡用牛刀。5. 避坑实录与验证方法5.1 溢出与取模大范围求和最常见的坑是中间结果溢出。n 10^12时n * (n 1) / 2大约是5 * 10^23远超 64 位整数范围必须每一步取模。即使题目要求输出精确值很多情况下也会因为结果过大而要求取模所以模板里最好默认模MOD只在明确需要精确值时去掉取模。另一个隐蔽的问题是r - l 1在整除分块中本身是O(√n)级别的数直接乘上S(n / l)也可能溢出正确的做法是把段长也取模后再乘。5.2 整除分块与取 id 的边界整除分块的r n / (n / l)是所有这类算法的地基。我见过不少人把n写成x或者在内层递归里忘了把n替换成当前函数的参数导致分块失效。Min_25 筛里get_id的映射更是高危区n / x的取整方向必须是向下取整id2的下标用的是n / x而不是n / w[t]写反一次就是访问越界。5.3 和对拍用线性筛验证小数据不管网上流传的模板多成熟新写的代码一定要对拍。我的做法是先写一个朴素版线性筛到10^6然后取n 10000、12345、99991等随机值跑同一组输入对比结果。对付μ和φ还可以直接用线性筛算出的前缀和做基准。对拍能筛掉九成以上的边界 bug尤其是 Min_25 筛递归边界和指数循环的错位。5.4 记忆化键与预处理阈值杜教筛的记忆化键尽量只用n本身因为递归中出现的所有参数都是⌊n/d⌋形式的整除块值数量有限。用unordered_map时提前reserve可以大幅减少重哈希。预处理阈值不是越大越好n^(2/3)只是经验值内存允许时稍微加大能减少递归深度但超过物理内存反而会因为换页拖慢整体速度。我一般以5 * 10^6到2 * 10^7为基准根据题目内存限制调节。最后的最后我在实际使用中的一个体会是杜教筛和 Min_25 筛虽然叫“筛法”但它们更像是“分治 整除分块 记忆化”的组合拳。你不需要把整套数论推导背得滚瓜烂熟但一定要亲手推一次 μ 和 φ 的递推式、写一次 Min_25 的第二段递归。第一次写完对拍通过之后你会觉得这类n 10^12的求和题不过如此。要是你只打算记住一个结论那就记住遇到大范围积性函数求和先看能不能构造卷积再看素数处的函数值是不是多项式前者走杜教筛后者走 Min_25 筛。如果你之后碰到自定义的奇怪积性函数也建议先用线性筛把1到10^5的数据打出来观察一下规律再决定用哪种筛法。很多“看似不可做”的题最后都倒在这两个模板面前。
