做算法竞赛题这几年数论里最让我头疼的其实不是快速幂、不是扩展欧几里得而是莫比乌斯函数和反演那一套。公式看得懂一到题目里就不知道往哪儿套。直到把某OJ上的第215题《破译密码》翻来覆去啃了几遍才真正把容斥、莫比乌斯函数和数论分块这三件事串成了一条线。这篇文章就把我这套完整的思考路径写出来不绕弯子从题目本身出发一步步推到你也能直接照着写的分块代码。适合刚学完数论基础、准备进阶莫比乌斯反演的选手也适合刷题时被这道题卡住、只知道抄模板但没想明白为什么的人。1. 题目本质把gcdd翻译成数论计数问题1.1 题意与数据规模题目描述是我见过的最简洁的那种给定三个正整数 a、b、d问有多少对整数 (x, y)满足 1 ≤ x ≤ a、1 ≤ y ≤ b并且 gcd(x, y) d。数据范围很关键a 和 b 最大到 50000查询次数 T 也最大到 50000。先别急着上手暴力我帮你算一笔账如果一次查询就双层 for 枚举 x 和 y50000 × 50000 2.5 × 10^9C 单次要跑好几秒而 T 可是有五万次。哪怕题目只给一次查询这个数据规模也足以劝退纯枚举做法。所以这题真正要解决的是怎么把单次查询从 O(ab) 降下来最终目标是 O(√min(a/d, b/d)) 这个量级。1.2 一个容易被忽略的特判d 必须是 a 和 b 的因子很多题解不会刻意提这点但我在实战中发现它很重要。题目要求 x 能被 d 整除所以 x 只能是 d 的倍数。在 1 到 a 之间d 的倍数一共有 a/d 向下取整个。如果 a 本身不是 d 的倍数比如 a 7、d 3那 a/d 2能选的 x 只有 3 和 6最大也就 6永远不会达到 a。这时如果还想着“用 a 除以 d 得到 n”n 其实只是整数部分的 2问题本身依然有解。但如果 a 比 d 小比如 a 2、d 5那 a/d 01 到 a 之间根本没有 d 的倍数答案直接就是 0。我当时就在这踩过坑提前把 a、b 都除以 d 得到 n、m然后发现 n 或 m 等于 0直接进入分块循环导致除零错误。正确的做法是先把 a 和 d、b 和 d 的关系看清楚a d 且 b d 只是基本条件更严谨的判断是如果 a d 或 b d或者更一般地a 与 b 对 d 整除后有一个为 0直接输出 0。后面的公式推导里还会看到n 0 会让分块循环无法正常执行。1.3 等价变换把 gcd(x, y) d 变成 gcd(x, y) 1既然 gcd(x, y) d那 x 一定是 d 的倍数y 也一定是 d 的倍数。把 x 写成 d·xy 写成 d·y原条件就等价变成了 gcd(d·x, d·y) d。两边约掉 d得到 gcd(x, y) 1也就是 x 和 y 互质。这一步没有任何技巧性纯粹是提取公因子但它把整个问题的目标重写了不再数“最大公约数恰好等于 d”的对数而是数“约分后互质”的对数。代价仅仅是让 x 的取值范围从 [1, a] 缩到 [1, n]其中 n a/dy 的范围变成 [1, m]其中 m b/d。如果 a 或 b 不是 d 的倍数那 n 或 m 其实是向下取整的结果不影响正确性因为能选的就是这么多。1.4 为什么互质计数依然不好直接做你可能会说互质不就是“没有公共质因子”吗那我枚举每个 x对每个 y 做一次 gcd 判断复杂度 O(nm)数据一大照样超时。而且互质的条件不像“被某个数整除”那样好数个数因为“没有公共质因子”是排除式条件排除对象是所有质数不是一个数。正是这个“排除”的特征天然指向容斥原理总数减去所有不互质的情况而每种不互质的情况都可以用“被某个约数同时整除”来刻画。2. 容斥视角下的莫比乌斯函数不背公式也能推出来2.1 用整除关系给集合计数设 F(t) 表示满足以下条件的 (x, y) 对数t 能同时整除 x 和 y也就是 t | gcd(x, y)。这个函数非常好算。在 1 到 n 的范围里t 的倍数有 n/t 向下取整个在 1 到 m 的范围里t 的倍数有 m/t 向下取整个。因为 x 和 y 的选择相互独立乘法原理直接得到F(t) ⌊n / t⌋ · ⌊m / t⌋比如 n m 4t 2 时x 可以选 2、4y 可以选 2、4所以 F(2) 2 × 2 4。2.2 容斥展开从总数里筛出互质对我们要的是 gcd(x, y) 1 的对数。总数显然是 F(1) n·m。但 F(1) 里包含了所有 gcd(x, y) 1 的对需要筛掉。怎么筛用容斥的经典套路。设所有质数依次为 p1, p2, ...。先把所有被某个质数 p 同时整除的对数减掉即减去 ∑ F(p)。但是如果一个数对同时能被 p 和 q 整除也就是 gcd 是 p 和 q 的倍数那它被减了两次要加回来即加回 ∑ F(p·q)。以此类推三个不同质数的乘积要减掉四个的乘积要加回来。于是答案可以展开成ans F(1) − ∑ F(p) ∑ F(p·q) − ∑ F(p·q·r) …注意展开式中出现的每一项对应的整数 t 都必须是一些互不相同的质数的乘积。如果 t 里面含平方因子例如 12 2² × 3它不会出现在这一串容斥里因为容斥只对“不同的质数集合”取交集同一个质数不会交两次。2.3 容斥系数就是莫比乌斯函数观察上面的展开式每一项 F(t) 前面的系数只取决于 t 的质因子构成t 1对应 F(1)系数为 1t 是 k 个互不相同质数的乘积系数为 (−1)^kt 含平方因子根本不出现等价于系数为 0。这就是莫比乌斯函数 μ(t) 的定义μ(t) 1当 t 1 μ(t) (−1)^k当 t 是 k 个互不相同质数的乘积 μ(t) 0当 t 含有平方因子。于是上面那一长串容斥可以压缩成一个求和ans ∑_{k1}^{min(n, m)} μ(k) · ⌊n / k⌋ · ⌊m / k⌋这里的求和上限取 min(n, m) 是因为当 k n 时 ⌊n/k⌋ 0当 k m 时 ⌊m/k⌋ 0F(k) 变成 0贡献为零。所以循环到 min(n, m) 就够了。2.4 这个公式和莫比乌斯反演的关系如果你学过尤拉筛和反演会看到教科书写着若 g(n) ∑_{d|n} f(d)则 f(n) ∑_{d|n} μ(d)·g(n/d)。这道题的公式其实可以看作反演的一种具体形式。但我觉得对做题来说反演公式更像是总结而容斥推导才是真正的理解入口。我自己更习惯这样记忆μ 函数是“容斥的系数表”它回答的问题只有一个——某个整数 t 在容斥里的加减号是什么。如果你在考场上忘了 μ 的具体定义临场画一个质因数分解看看有几个不同的质因子再检查有没有平方因子也能手算出 μ 值。这个习惯帮我避免了很多次“公式背反了”的尴尬。3. 莫比乌斯函数的线性筛边界细节决定成败3.1 为什么必须用筛法预处理从公式看一次查询需要枚举 k 从 1 到 min(n, m)最多 50000。如果每个 k 都现场质因数分解来算 μ(k)一次查询的成本就会带上 log 因子五万次查询照样会卡。正确的做法是提前用线性筛预处理 1 到 50000 的所有 μ 值顺便算出前缀和。前缀和是为了给后面的数论分块用这个下一节会讲。3.2 线性筛求 μ 的完整模板线性筛能求 μ 的核心在于每个合数只被它的最小质因子筛到而且筛的过程中可以顺便得知新数的质因子构成。模板代码我贴在下面这是我在多个题解基础上精简过的版本#include bits/stdc.h using namespace std; typedef long long ll; const int MAXN 50000 5; int primes[MAXN], cnt; bool vis[MAXN]; ll mu[MAXN]; ll sumMu[MAXN]; void sieve(int n) { mu[1] 1; for (int i 2; i n; i) { if (!vis[i]) { primes[cnt] i; mu[i] -1; } for (int j 1; j cnt i * primes[j] n; j) { vis[i * primes[j]] true; if (i % primes[j] 0) { mu[i * primes[j]] 0; break; } else { mu[i * primes[j]] -mu[i]; } } } for (int i 1; i n; i) { sumMu[i] sumMu[i - 1] mu[i]; } }这段代码里最需要注意的地方有三个第一mu[1] 必须单独赋值为 1因为循环从 i 2 开始。这一点漏了后面所有前缀和都会偏。第二当 i % primes[j] 0 时说明 i 已经包含质因子 primes[j]那么新的数 i × primes[j] 里primes[j] 的指数至少为 2也就是含平方因子所以 μ 直接置 0然后立刻 break。这个 break 不仅是为了正确性更是线性筛“每个数只被筛一次”的关键。第三在 i % primes[j] ! 0 的分支里i 不含 primes[j]所以 i × primes[j] 相当于在 i 的基础上新增了一个不同质因子质因子数量加 1符号自然取反也就是 μ[i × primes[j]] -μ[i]。3.3 手算几个值验证筛法我平时写完筛子都会随手验几个数μ(1) 1μ(2) -1μ(3) -1μ(4) 0μ(5) -1μ(6) 1μ(8) 0μ(30) (-1)^3 -1。如果筛法输出和这些一致基本可以认为代码正确。3.4 一个值得注意的性能点MAXN 不要开得刚刚好。题目给 a、b 最大 50000那就开 50005 或者更大一点的常量。如果你在练更进阶的版本数据范围到 1e7注意 int 乘法的溢出i * primes[j] 可能超过 int 上界要转成 long long 判断。我上面的模板里就是这样处理的虽然 50000 这个范围不会溢出但养成习惯总没坏处。4. 数论分块优化为什么答案能O(√n)计算4.1 观察商数的成段性质我们现在要算的是ans ∑_{k1}^{min(n,m)} μ(k) · ⌊n / k⌋ · ⌊m / k⌋如果直接 for 循环i 从 1 到 min(n,m)单次查询复杂度是 O(min(n,m))。看起来好像也能接受50000 和一次查询还行但有 50000 次查询直接变成 2.5 × 10^9 次迭代还是不行。关键观察⌊n/k⌋ 的值并不会随着 k 每个 1 都变化。比如 n 10k: 1 2 3 4 5 6 7 8 9 10 ⌊10/k⌋: 10 5 3 3 2 2 1 1 1 1可以发现值是一段一段相同的。而且这样的段数很少大约只有 2√n 个。这就是数论分块也叫整除分块能提速的根本原因。4.2 怎么快速求一段区间的右端点假设当前区间左端点是 l令 q ⌊n/l⌋我们希望找到最大的 r使得 ⌊n/k⌋ q 对所有 k ∈ [l, r] 都成立。答案是r ⌊n / q⌋为什么因为当 k ⌊n/q⌋ 时⌊n/k⌋ 会小于 q。这个结论可以直接背下来也可以这样理解⌊n/k⌋ ≥ q 等价于 k ≤ ⌊n/q⌋而恰好 ⌊n/k⌋ q 的范围就是 [l, ⌊n/q⌋]。这是整除分块最核心的公式。4.3 双变量分块右端点取最小值本题有两个变量 n 和 m。我们希望在一个区间 [l, r] 内⌊n/k⌋ 和 ⌊m/k⌋ 都保持不变这样 ⌊n/k⌋·⌊m/k⌋ 作为一个整体也不用变。做法是分别求右端点r1 n / (n / l)r2 m / (m / l)然后取 r min(r1, r2)。这样在 [l, r] 内两个商数都不变我们就可以安心把这段区间上的 μ(k) 加起来然后一次性乘上不变的两个商。4.4 区间 μ 求和改成前缀和既然一段区间里的两个商数都一样那这一整段的贡献就是( sumMu[r] - sumMu[l-1] ) × ⌊n/l⌋ × ⌊m/l⌋这里 sumMu 就是上一节线性筛算出的莫比乌斯函数前缀和。整个算法的跳跃次数等于 n 的分块段数和 m 的分块段数中较小的数量级大约是 2√min(n,m)。分块代码ll solve(int n, int m) { ll ans 0; int limit min(n, m); for (int l 1, r; l limit; l r 1) { int qn n / l; int qm m / l; r min(n / qn, m / qm); ans (sumMu[r] - sumMu[l - 1]) * qn * qm; } return ans; }注意循环条件里必须写 l limit而不是 l n。因为当 l 超过 m 之后qm m / l 会变成 0那么 m / qm 就是除零错误。这个坑我在比赛里真踩过当时还以为是数据有问题后来才发现是自己的循环边界写错了。4.5 分块复杂度实测我拿极限数据 a b 50000、d 1、T 50000 实测过预处理一次筛法之后solve 函数每个查询的分块段数大约只有 440 个左右整体运行时间在几百毫秒内。相比 O(min(n,m)) 的 2.5 × 10^9 次迭代这个优化量级完全不在一个档次。所以“数论分块”这道题的正解不是可有可无的优化而是必须做的核心步骤。5. 完整代码与边界情况复盘5.1 完整可提交代码把筛法、分块查询、输入输出整合起来完整代码如下#include bits/stdc.h using namespace std; typedef long long ll; const int MAXN 50000 5; int primes[MAXN], cnt; bool vis[MAXN]; ll mu[MAXN]; ll sumMu[MAXN]; void initMobius(int n) { mu[1] 1; for (int i 2; i n; i) { if (!vis[i]) { primes[cnt] i; mu[i] -1; } for (int j 1; j cnt (ll)i * primes[j] n; j) { vis[i * primes[j]] true; if (i % primes[j] 0) { mu[i * primes[j]] 0; break; } else { mu[i * primes[j]] -mu[i]; } } } for (int i 1; i n; i) { sumMu[i] sumMu[i - 1] mu[i]; } } ll calc(int n, int m) { ll ans 0; int limit min(n, m); for (int l 1, r; l limit; l r 1) { int qn n / l; int qm m / l; r min(n / qn, m / qm); ans (sumMu[r] - sumMu[l - 1]) * (ll)qn * (ll)qm; } return ans; } int main() { ios::sync_with_stdio(false); cin.tie(nullptr); initMobius(MAXN - 1); int T; cin T; while (T--) { int a, b, d; cin a b d; if (a d || b d) { cout 0 \n; continue; } int n a / d; int m b / d; cout calc(n, m) \n; } return 0; }这个版本里sumMu 和 mu 用 long long 是因为后面乘法可能出现 qn × qm 最大约 2.5 × 10^9如果再用 int 中途截断再转 long long会在赋值前就已经溢出。直接全链用 long long 最省心。5.2 手算一遍验证正确性用题目最简单的例子a 4, b 4, d 1也就是 n m 4。根据公式展开k 1μ(1) 1⌊4/1⌋·⌊4/1⌋ 16贡献 16 k 2μ(2) -1⌊4/2⌋·⌊4/2⌋ 4贡献 -4 k 3μ(3) -1⌊4/3⌋·⌊4/3⌋ 1贡献 -1 k 4μ(4) 0贡献 0加起来 ans 16 - 4 - 1 11。实际枚举 16 对互质对是 (1,1)、(1,2)、(1,3)、(1,4)、(2,1)、(2,3)、(3,1)、(3,2)、(3,4)、(4,1)、(4,3)正好 11 对。剩下的 5 对 (2,2)、(2,4)、(3,3)、(4,2)、(4,4)gcd 都大于 1。公式和枚举对得上。5.3 多组数据的预处理只需要一次initMobius 放在所有查询之前执行一次而不是每个查询都跑。这一点很基础但容易在紧张时写错。筛法复杂度 O(MAXN)即便 MAXN 开到 50000也只需要准备一次后面每个查询都是 O(√min(n,m))。6. 实战踩坑取整陷阱与复杂度误区6.1 前缀和数组类型引发的错误我第一次提交时把 sumMu 定义成 int结果在大范围数据下 WA 了几发。原因在于 sumMu 的值本身最大也就三万多int 存得下。但计算贡献时乘法 (sumMu[r] - sumMu[l-1]) * qn * qm 在某个区间 qn50000、qm50000、sumMu 差值为 1 的情况下结果能达到 2.5 × 10^9已经超出 int 范围。编译器按 int 算完再赋给 long long溢出已经发生。所以要么在乘法前显式把每个操作数转成 long long要么直接把 sumMu、mu 定义成 long long。后一种更省事。6.2 分块右端点公式写错正确公式是 r min(n / (n / l), m / (m / l))。我见过初学者写成 r min(n, m) / (min(n, m) / l)这完全不对。因为两个分式的除数不同不能简单合并。还有一种错误是把 r1 写成 n / (n / l) 后就不再管 m直接让区间按 n 的段走。这样如果 m 在中间某个点先发生跳变那 [l, r1] 这段区间里 qm 已经变了你还按之前的 qm 乘结果必错。调试这种错误很痛苦因为它只在特定数据下出错。6.3 循环边界忘记取 min导致除零循环条件写成 for (l 1; l n; l r1) 是致命的。当 l 超过 m 之后qm m / l 0然后 m / qm 就除零了。而且一旦 l 越过 m后续所有右端点都会出问题。我在第一次独立写双变量分块时就犯了这个错排查了很久才发现最后加了一个 limit min(n, m) 才解决。以后凡是双变量分块我都在代码注释里标注“边界必须 min”。6.4 忘记特判 a d 或 b d如果不特判n 或 m 可能为 0for 循环一开始 limit 0直接跳过循环答案输出 0。看起来也能过但这不是理所当然的。更隐性的情况是a 可以是 d 的倍数但 b 不是这时 m b/d 是向下取整问题本身仍然有解不能输出 0。所以要分清什么时候是真的无解a d 或 b d什么时候只是上界缩小。6.5 复杂度误区O(√n) 的区间数量究竟是多少很多人误以为分块后复杂度就是 O(√n)其实准确说是 O(2√n) 量级。对于单变量 n当 k ≤ √n 时⌊n/k⌋ 有大约 √n 个取值当 k √n 时⌊n/k⌋ 的取值范围在 1 到 √n 之间也是大约 √n 个取值。两部分加起来不超过 2√n 段。双变量分块的段数不超过两个单变量段数之和。所以整体复杂度可以放心写 O(√min(n,m)) 或 O(√n √m)。最后想说的一点体会莫比乌斯函数这玩意刚接触时总觉得像天书符号乱、公式多。但《破译密码》这道题给了我一个非常好的锚点它把容斥、μ 函数、分块三件事揉在一起却又没有故意设置障碍。我把推导过程完完整整写下来之后最大的收获是认识到“μ 值就是容斥系数”这句话远比“莫比乌斯反演公式”实用。遇到计数题先想想能不能用容斥展开再想想展开系数是不是正好落在 μ 上最后看能不能分块加速这个思考顺序比硬背公式可靠得多。如果你也被这道题卡过希望这篇记录能让你少走一点弯路。
做算法竞赛题这几年数论里最让我头疼的其实不是快速幂、不是扩展欧几里得而是莫比乌斯函数和反演那一套。公式看得懂一到题目里就不知道往哪儿套。直到把某OJ上的第215题《破译密码》翻来覆去啃了几遍才真正把容斥、莫比乌斯函数和数论分块这三件事串成了一条线。这篇文章就把我这套完整的思考路径写出来不绕弯子从题目本身出发一步步推到你也能直接照着写的分块代码。适合刚学完数论基础、准备进阶莫比乌斯反演的选手也适合刷题时被这道题卡住、只知道抄模板但没想明白为什么的人。1. 题目本质把gcdd翻译成数论计数问题1.1 题意与数据规模题目描述是我见过的最简洁的那种给定三个正整数 a、b、d问有多少对整数 (x, y)满足 1 ≤ x ≤ a、1 ≤ y ≤ b并且 gcd(x, y) d。数据范围很关键a 和 b 最大到 50000查询次数 T 也最大到 50000。先别急着上手暴力我帮你算一笔账如果一次查询就双层 for 枚举 x 和 y50000 × 50000 2.5 × 10^9C 单次要跑好几秒而 T 可是有五万次。哪怕题目只给一次查询这个数据规模也足以劝退纯枚举做法。所以这题真正要解决的是怎么把单次查询从 O(ab) 降下来最终目标是 O(√min(a/d, b/d)) 这个量级。1.2 一个容易被忽略的特判d 必须是 a 和 b 的因子很多题解不会刻意提这点但我在实战中发现它很重要。题目要求 x 能被 d 整除所以 x 只能是 d 的倍数。在 1 到 a 之间d 的倍数一共有 a/d 向下取整个。如果 a 本身不是 d 的倍数比如 a 7、d 3那 a/d 2能选的 x 只有 3 和 6最大也就 6永远不会达到 a。这时如果还想着“用 a 除以 d 得到 n”n 其实只是整数部分的 2问题本身依然有解。但如果 a 比 d 小比如 a 2、d 5那 a/d 01 到 a 之间根本没有 d 的倍数答案直接就是 0。我当时就在这踩过坑提前把 a、b 都除以 d 得到 n、m然后发现 n 或 m 等于 0直接进入分块循环导致除零错误。正确的做法是先把 a 和 d、b 和 d 的关系看清楚a d 且 b d 只是基本条件更严谨的判断是如果 a d 或 b d或者更一般地a 与 b 对 d 整除后有一个为 0直接输出 0。后面的公式推导里还会看到n 0 会让分块循环无法正常执行。1.3 等价变换把 gcd(x, y) d 变成 gcd(x, y) 1既然 gcd(x, y) d那 x 一定是 d 的倍数y 也一定是 d 的倍数。把 x 写成 d·xy 写成 d·y原条件就等价变成了 gcd(d·x, d·y) d。两边约掉 d得到 gcd(x, y) 1也就是 x 和 y 互质。这一步没有任何技巧性纯粹是提取公因子但它把整个问题的目标重写了不再数“最大公约数恰好等于 d”的对数而是数“约分后互质”的对数。代价仅仅是让 x 的取值范围从 [1, a] 缩到 [1, n]其中 n a/dy 的范围变成 [1, m]其中 m b/d。如果 a 或 b 不是 d 的倍数那 n 或 m 其实是向下取整的结果不影响正确性因为能选的就是这么多。1.4 为什么互质计数依然不好直接做你可能会说互质不就是“没有公共质因子”吗那我枚举每个 x对每个 y 做一次 gcd 判断复杂度 O(nm)数据一大照样超时。而且互质的条件不像“被某个数整除”那样好数个数因为“没有公共质因子”是排除式条件排除对象是所有质数不是一个数。正是这个“排除”的特征天然指向容斥原理总数减去所有不互质的情况而每种不互质的情况都可以用“被某个约数同时整除”来刻画。2. 容斥视角下的莫比乌斯函数不背公式也能推出来2.1 用整除关系给集合计数设 F(t) 表示满足以下条件的 (x, y) 对数t 能同时整除 x 和 y也就是 t | gcd(x, y)。这个函数非常好算。在 1 到 n 的范围里t 的倍数有 n/t 向下取整个在 1 到 m 的范围里t 的倍数有 m/t 向下取整个。因为 x 和 y 的选择相互独立乘法原理直接得到F(t) ⌊n / t⌋ · ⌊m / t⌋比如 n m 4t 2 时x 可以选 2、4y 可以选 2、4所以 F(2) 2 × 2 4。2.2 容斥展开从总数里筛出互质对我们要的是 gcd(x, y) 1 的对数。总数显然是 F(1) n·m。但 F(1) 里包含了所有 gcd(x, y) 1 的对需要筛掉。怎么筛用容斥的经典套路。设所有质数依次为 p1, p2, ...。先把所有被某个质数 p 同时整除的对数减掉即减去 ∑ F(p)。但是如果一个数对同时能被 p 和 q 整除也就是 gcd 是 p 和 q 的倍数那它被减了两次要加回来即加回 ∑ F(p·q)。以此类推三个不同质数的乘积要减掉四个的乘积要加回来。于是答案可以展开成ans F(1) − ∑ F(p) ∑ F(p·q) − ∑ F(p·q·r) …注意展开式中出现的每一项对应的整数 t 都必须是一些互不相同的质数的乘积。如果 t 里面含平方因子例如 12 2² × 3它不会出现在这一串容斥里因为容斥只对“不同的质数集合”取交集同一个质数不会交两次。2.3 容斥系数就是莫比乌斯函数观察上面的展开式每一项 F(t) 前面的系数只取决于 t 的质因子构成t 1对应 F(1)系数为 1t 是 k 个互不相同质数的乘积系数为 (−1)^kt 含平方因子根本不出现等价于系数为 0。这就是莫比乌斯函数 μ(t) 的定义μ(t) 1当 t 1 μ(t) (−1)^k当 t 是 k 个互不相同质数的乘积 μ(t) 0当 t 含有平方因子。于是上面那一长串容斥可以压缩成一个求和ans ∑_{k1}^{min(n, m)} μ(k) · ⌊n / k⌋ · ⌊m / k⌋这里的求和上限取 min(n, m) 是因为当 k n 时 ⌊n/k⌋ 0当 k m 时 ⌊m/k⌋ 0F(k) 变成 0贡献为零。所以循环到 min(n, m) 就够了。2.4 这个公式和莫比乌斯反演的关系如果你学过尤拉筛和反演会看到教科书写着若 g(n) ∑_{d|n} f(d)则 f(n) ∑_{d|n} μ(d)·g(n/d)。这道题的公式其实可以看作反演的一种具体形式。但我觉得对做题来说反演公式更像是总结而容斥推导才是真正的理解入口。我自己更习惯这样记忆μ 函数是“容斥的系数表”它回答的问题只有一个——某个整数 t 在容斥里的加减号是什么。如果你在考场上忘了 μ 的具体定义临场画一个质因数分解看看有几个不同的质因子再检查有没有平方因子也能手算出 μ 值。这个习惯帮我避免了很多次“公式背反了”的尴尬。3. 莫比乌斯函数的线性筛边界细节决定成败3.1 为什么必须用筛法预处理从公式看一次查询需要枚举 k 从 1 到 min(n, m)最多 50000。如果每个 k 都现场质因数分解来算 μ(k)一次查询的成本就会带上 log 因子五万次查询照样会卡。正确的做法是提前用线性筛预处理 1 到 50000 的所有 μ 值顺便算出前缀和。前缀和是为了给后面的数论分块用这个下一节会讲。3.2 线性筛求 μ 的完整模板线性筛能求 μ 的核心在于每个合数只被它的最小质因子筛到而且筛的过程中可以顺便得知新数的质因子构成。模板代码我贴在下面这是我在多个题解基础上精简过的版本#include bits/stdc.h using namespace std; typedef long long ll; const int MAXN 50000 5; int primes[MAXN], cnt; bool vis[MAXN]; ll mu[MAXN]; ll sumMu[MAXN]; void sieve(int n) { mu[1] 1; for (int i 2; i n; i) { if (!vis[i]) { primes[cnt] i; mu[i] -1; } for (int j 1; j cnt i * primes[j] n; j) { vis[i * primes[j]] true; if (i % primes[j] 0) { mu[i * primes[j]] 0; break; } else { mu[i * primes[j]] -mu[i]; } } } for (int i 1; i n; i) { sumMu[i] sumMu[i - 1] mu[i]; } }这段代码里最需要注意的地方有三个第一mu[1] 必须单独赋值为 1因为循环从 i 2 开始。这一点漏了后面所有前缀和都会偏。第二当 i % primes[j] 0 时说明 i 已经包含质因子 primes[j]那么新的数 i × primes[j] 里primes[j] 的指数至少为 2也就是含平方因子所以 μ 直接置 0然后立刻 break。这个 break 不仅是为了正确性更是线性筛“每个数只被筛一次”的关键。第三在 i % primes[j] ! 0 的分支里i 不含 primes[j]所以 i × primes[j] 相当于在 i 的基础上新增了一个不同质因子质因子数量加 1符号自然取反也就是 μ[i × primes[j]] -μ[i]。3.3 手算几个值验证筛法我平时写完筛子都会随手验几个数μ(1) 1μ(2) -1μ(3) -1μ(4) 0μ(5) -1μ(6) 1μ(8) 0μ(30) (-1)^3 -1。如果筛法输出和这些一致基本可以认为代码正确。3.4 一个值得注意的性能点MAXN 不要开得刚刚好。题目给 a、b 最大 50000那就开 50005 或者更大一点的常量。如果你在练更进阶的版本数据范围到 1e7注意 int 乘法的溢出i * primes[j] 可能超过 int 上界要转成 long long 判断。我上面的模板里就是这样处理的虽然 50000 这个范围不会溢出但养成习惯总没坏处。4. 数论分块优化为什么答案能O(√n)计算4.1 观察商数的成段性质我们现在要算的是ans ∑_{k1}^{min(n,m)} μ(k) · ⌊n / k⌋ · ⌊m / k⌋如果直接 for 循环i 从 1 到 min(n,m)单次查询复杂度是 O(min(n,m))。看起来好像也能接受50000 和一次查询还行但有 50000 次查询直接变成 2.5 × 10^9 次迭代还是不行。关键观察⌊n/k⌋ 的值并不会随着 k 每个 1 都变化。比如 n 10k: 1 2 3 4 5 6 7 8 9 10 ⌊10/k⌋: 10 5 3 3 2 2 1 1 1 1可以发现值是一段一段相同的。而且这样的段数很少大约只有 2√n 个。这就是数论分块也叫整除分块能提速的根本原因。4.2 怎么快速求一段区间的右端点假设当前区间左端点是 l令 q ⌊n/l⌋我们希望找到最大的 r使得 ⌊n/k⌋ q 对所有 k ∈ [l, r] 都成立。答案是r ⌊n / q⌋为什么因为当 k ⌊n/q⌋ 时⌊n/k⌋ 会小于 q。这个结论可以直接背下来也可以这样理解⌊n/k⌋ ≥ q 等价于 k ≤ ⌊n/q⌋而恰好 ⌊n/k⌋ q 的范围就是 [l, ⌊n/q⌋]。这是整除分块最核心的公式。4.3 双变量分块右端点取最小值本题有两个变量 n 和 m。我们希望在一个区间 [l, r] 内⌊n/k⌋ 和 ⌊m/k⌋ 都保持不变这样 ⌊n/k⌋·⌊m/k⌋ 作为一个整体也不用变。做法是分别求右端点r1 n / (n / l)r2 m / (m / l)然后取 r min(r1, r2)。这样在 [l, r] 内两个商数都不变我们就可以安心把这段区间上的 μ(k) 加起来然后一次性乘上不变的两个商。4.4 区间 μ 求和改成前缀和既然一段区间里的两个商数都一样那这一整段的贡献就是( sumMu[r] - sumMu[l-1] ) × ⌊n/l⌋ × ⌊m/l⌋这里 sumMu 就是上一节线性筛算出的莫比乌斯函数前缀和。整个算法的跳跃次数等于 n 的分块段数和 m 的分块段数中较小的数量级大约是 2√min(n,m)。分块代码ll solve(int n, int m) { ll ans 0; int limit min(n, m); for (int l 1, r; l limit; l r 1) { int qn n / l; int qm m / l; r min(n / qn, m / qm); ans (sumMu[r] - sumMu[l - 1]) * qn * qm; } return ans; }注意循环条件里必须写 l limit而不是 l n。因为当 l 超过 m 之后qm m / l 会变成 0那么 m / qm 就是除零错误。这个坑我在比赛里真踩过当时还以为是数据有问题后来才发现是自己的循环边界写错了。4.5 分块复杂度实测我拿极限数据 a b 50000、d 1、T 50000 实测过预处理一次筛法之后solve 函数每个查询的分块段数大约只有 440 个左右整体运行时间在几百毫秒内。相比 O(min(n,m)) 的 2.5 × 10^9 次迭代这个优化量级完全不在一个档次。所以“数论分块”这道题的正解不是可有可无的优化而是必须做的核心步骤。5. 完整代码与边界情况复盘5.1 完整可提交代码把筛法、分块查询、输入输出整合起来完整代码如下#include bits/stdc.h using namespace std; typedef long long ll; const int MAXN 50000 5; int primes[MAXN], cnt; bool vis[MAXN]; ll mu[MAXN]; ll sumMu[MAXN]; void initMobius(int n) { mu[1] 1; for (int i 2; i n; i) { if (!vis[i]) { primes[cnt] i; mu[i] -1; } for (int j 1; j cnt (ll)i * primes[j] n; j) { vis[i * primes[j]] true; if (i % primes[j] 0) { mu[i * primes[j]] 0; break; } else { mu[i * primes[j]] -mu[i]; } } } for (int i 1; i n; i) { sumMu[i] sumMu[i - 1] mu[i]; } } ll calc(int n, int m) { ll ans 0; int limit min(n, m); for (int l 1, r; l limit; l r 1) { int qn n / l; int qm m / l; r min(n / qn, m / qm); ans (sumMu[r] - sumMu[l - 1]) * (ll)qn * (ll)qm; } return ans; } int main() { ios::sync_with_stdio(false); cin.tie(nullptr); initMobius(MAXN - 1); int T; cin T; while (T--) { int a, b, d; cin a b d; if (a d || b d) { cout 0 \n; continue; } int n a / d; int m b / d; cout calc(n, m) \n; } return 0; }这个版本里sumMu 和 mu 用 long long 是因为后面乘法可能出现 qn × qm 最大约 2.5 × 10^9如果再用 int 中途截断再转 long long会在赋值前就已经溢出。直接全链用 long long 最省心。5.2 手算一遍验证正确性用题目最简单的例子a 4, b 4, d 1也就是 n m 4。根据公式展开k 1μ(1) 1⌊4/1⌋·⌊4/1⌋ 16贡献 16 k 2μ(2) -1⌊4/2⌋·⌊4/2⌋ 4贡献 -4 k 3μ(3) -1⌊4/3⌋·⌊4/3⌋ 1贡献 -1 k 4μ(4) 0贡献 0加起来 ans 16 - 4 - 1 11。实际枚举 16 对互质对是 (1,1)、(1,2)、(1,3)、(1,4)、(2,1)、(2,3)、(3,1)、(3,2)、(3,4)、(4,1)、(4,3)正好 11 对。剩下的 5 对 (2,2)、(2,4)、(3,3)、(4,2)、(4,4)gcd 都大于 1。公式和枚举对得上。5.3 多组数据的预处理只需要一次initMobius 放在所有查询之前执行一次而不是每个查询都跑。这一点很基础但容易在紧张时写错。筛法复杂度 O(MAXN)即便 MAXN 开到 50000也只需要准备一次后面每个查询都是 O(√min(n,m))。6. 实战踩坑取整陷阱与复杂度误区6.1 前缀和数组类型引发的错误我第一次提交时把 sumMu 定义成 int结果在大范围数据下 WA 了几发。原因在于 sumMu 的值本身最大也就三万多int 存得下。但计算贡献时乘法 (sumMu[r] - sumMu[l-1]) * qn * qm 在某个区间 qn50000、qm50000、sumMu 差值为 1 的情况下结果能达到 2.5 × 10^9已经超出 int 范围。编译器按 int 算完再赋给 long long溢出已经发生。所以要么在乘法前显式把每个操作数转成 long long要么直接把 sumMu、mu 定义成 long long。后一种更省事。6.2 分块右端点公式写错正确公式是 r min(n / (n / l), m / (m / l))。我见过初学者写成 r min(n, m) / (min(n, m) / l)这完全不对。因为两个分式的除数不同不能简单合并。还有一种错误是把 r1 写成 n / (n / l) 后就不再管 m直接让区间按 n 的段走。这样如果 m 在中间某个点先发生跳变那 [l, r1] 这段区间里 qm 已经变了你还按之前的 qm 乘结果必错。调试这种错误很痛苦因为它只在特定数据下出错。6.3 循环边界忘记取 min导致除零循环条件写成 for (l 1; l n; l r1) 是致命的。当 l 超过 m 之后qm m / l 0然后 m / qm 就除零了。而且一旦 l 越过 m后续所有右端点都会出问题。我在第一次独立写双变量分块时就犯了这个错排查了很久才发现最后加了一个 limit min(n, m) 才解决。以后凡是双变量分块我都在代码注释里标注“边界必须 min”。6.4 忘记特判 a d 或 b d如果不特判n 或 m 可能为 0for 循环一开始 limit 0直接跳过循环答案输出 0。看起来也能过但这不是理所当然的。更隐性的情况是a 可以是 d 的倍数但 b 不是这时 m b/d 是向下取整问题本身仍然有解不能输出 0。所以要分清什么时候是真的无解a d 或 b d什么时候只是上界缩小。6.5 复杂度误区O(√n) 的区间数量究竟是多少很多人误以为分块后复杂度就是 O(√n)其实准确说是 O(2√n) 量级。对于单变量 n当 k ≤ √n 时⌊n/k⌋ 有大约 √n 个取值当 k √n 时⌊n/k⌋ 的取值范围在 1 到 √n 之间也是大约 √n 个取值。两部分加起来不超过 2√n 段。双变量分块的段数不超过两个单变量段数之和。所以整体复杂度可以放心写 O(√min(n,m)) 或 O(√n √m)。最后想说的一点体会莫比乌斯函数这玩意刚接触时总觉得像天书符号乱、公式多。但《破译密码》这道题给了我一个非常好的锚点它把容斥、μ 函数、分块三件事揉在一起却又没有故意设置障碍。我把推导过程完完整整写下来之后最大的收获是认识到“μ 值就是容斥系数”这句话远比“莫比乌斯反演公式”实用。遇到计数题先想想能不能用容斥展开再想想展开系数是不是正好落在 μ 上最后看能不能分块加速这个思考顺序比硬背公式可靠得多。如果你也被这道题卡过希望这篇记录能让你少走一点弯路。
