从最小公倍数与容斥原理求解Codeforces Strange Function数论问题

1. 问题引入当“奇怪函数”遇上数论竞赛最近在整理一些经典的编程竞赛题目时又翻到了这道 Codeforces Round #743 (Div. 2) 的 C 题 “Strange Function”。题目本身描述很简洁定义一个函数 f(n)它返回最小的正整数 x使得 x 不是 n 的因子。例如f(1)2因为 1 能被 1 整除但 2 不能被 1 整除f(2)3因为 1 和 2 都是 2 的因子但 3 不是f(4)3因为 1 和 2 是 4 的因子但 3 不是。题目要求计算的是 S(n) f(1) f(2) ... f(n) 对某个大质数通常是 1e97取模的结果。n 的范围可以非常大最大到 1e16。看到这个数据范围任何试图暴力计算 f(n) 再累加的 O(n) 想法都可以直接放弃了。这明摆着是一道需要挖掘数学性质寻找高效计算方法的数论题。我第一次看到这个函数定义时也觉得它有点“奇怪”。但竞赛题的精妙之处就在于它往往从一个看似简单甚至有些“怪”的定义出发引导你深入思考数论中的核心概念。这道题的核心恰恰落在了最小公倍数 (LCM)和容斥原理这两个经典工具上。很多选手看到“容斥”可能会想到集合计数但在这里它以一种非常巧妙的方式帮助我们计算了在巨大范围内满足特定整除性质的数的个数。接下来我们就一步步拆解看看如何将天马行空的求和问题转化为可高效计算的数学模型。2. 函数 f(n) 的本质分析与转化思路直接计算每个 f(n) 是行不通的我们必须找到 f(n) 的另一种表达方式或者找到 S(n) 的整体计算公式。让我们先深入观察 f(n) 这个函数。根据定义f(n) 是最小的正整数 x使得 x 不是 n 的因子。换句话说对于所有比 x 小的正整数 k (1 ≤ k x)都有 k 是 n 的因子。这意味着什么这意味着 n 必须能被 1, 2, 3, ..., x-1 所有这些数同时整除。一个数能被一系列数同时整除等价于它能被这些数的最小公倍数 (LCM)整除。所以条件 “1, 2, ..., x-1 都是 n 的因子” 等价于 “LCM(1, 2, ..., x-1) 是 n 的因子”。而 f(n) x 意味着LCM(1,2,...,x-1) 能整除 n但 LCM(1,2,...,x) 不能整除 n。这给了我们一个全新的视角f(n) 的值取决于 n 能被前多少个连续正整数的 LCM 整除。更具体地说如果 n 能被 LCM(1) 1 整除显然总是成立那么 f(n) 至少为 2。如果 n 还能被 LCM(1,2)2 整除那么 f(n) 至少为 3。如果 n 还能被 LCM(1,2,3)6 整除那么 f(n) 至少为 4。以此类推。假设我们定义 L(k) LCM(1, 2, ..., k)。那么对于给定的 nf(n) 就是满足L(x-1) | n且L(x) ∤ n的那个 x。符号|表示整除∤表示不整除。因此计算 f(n) 的过程变成了检查 n 能被 L(1), L(2), L(3), ... 整除到哪一级。那么对于求和 S(n) Σ f(i)我们可以换个角度思考有多少个 i (1 ≤ i ≤ n)其 f(i) 的值等于某个特定的 x如果 f(i) x根据上述分析这意味着i 能被 L(x-1) 整除。i 不能被 L(x) 整除。满足条件1的 i 的个数是floor(n / L(x-1))即在 1 到 n 中 L(x-1) 的倍数的个数。 但是我们需要从中剔除那些同时也被 L(x) 整除的数因为那些数的 f(i) 至少是 x1。 而同时被 L(x-1) 和 L(x) 整除的数其实就是被 LCM(L(x-1), L(x)) 整除的数。由于 L(x) LCM(L(x-1), x)且 L(x-1) 显然整除 L(x)所以这个 LCM 就是 L(x) 本身。 因此满足f(i) x的 i 的个数就等于能被 L(x-1) 整除但不能被 L(x) 整除的数的个数。这正好是容斥原理的经典应用场景集合 A {i | L(x-1) | i}集合 B {i | L(x) | i}。我们要求的是在 A 中但不在 B 中的元素个数即 |A| - |A ∩ B|。由于 B 是 A 的子集L(x)是L(x-1)的倍数所以 |A ∩ B| |B|。 因此count(x) floor(n / L(x-1)) - floor(n / L(x))。于是总和 S(n) 可以表示为S(n) Σ_{x2}^{∞} [ x * count(x) ] Σ_{x2}^{∞} [ x * (floor(n / L(x-1)) - floor(n / L(x))) ]这里 x 的上限是无穷但实际上当 L(x-1) n 时floor(n / L(x-1)) 0后续项都为 0所以求和是有限的。我们只需要计算到 L(x-1) ≤ n 的那些 x 即可。3. 关键数学工具L(k) 序列的性质与高效生成上面的推导将问题转化为计算 L(k) LCM(1, 2, ..., k) 的序列以及对应的求和。L(k) 有一个非常重要的性质它增长得非常快。L(k) 是 1 到 k 所有数的最小公倍数它等于所有质数 p ≤ k 的 p 的最高次幂的乘积。例如L(1) 1L(2) LCM(1,2) 2L(3) LCM(1,2,3) 6L(4) LCM(1,2,3,4) LCM(6, 4) 12 因为 4 2^2而 L(3)62*3所以需要升级2的幂次L(5) LCM(12, 5) 60L(6) LCM(60, 6) 60 因为623已包含在602^23*5中L(7) 420L(8) LCM(420, 8) LCM(420, 2^3) 840 2的幂次从2升级到3观察可知L(k) 并非每次 k 增加1都会变化。只有当新的 k 引入了新的质因子或者提高了某个现有质因子的幂次时L(k) 才会更新。具体来说如果 k 是质数那么 L(k) L(k-1) * k。如果 k 是质数的幂次比如 p^e且 L(k-1) 中 p 的幂次小于 e那么 L(k) L(k-1) * p。否则L(k) L(k-1)。由于 n 最大为 1e16我们需要知道 L(k) 大概在 k 为多少时会超过这个值。通过计算或查阅已知数列OEIS A003418我们可以发现 L(1)1, L(2)2, L(3)6, L(4)12, L(5)60, L(6)60, L(7)420, L(8)840, L(9)2520, L(10)2520, L(11)27720, L(12)27720, L(13)360360, L(14)360360, L(15)360360, L(16)720720, L(17)12252240, L(18)12252240, L(19)232792560, L(20)232792560, L(21)232792560, L(22)232792560, L(23)5354228880, L(24)5354228880, L(25)26771144400, L(26)26771144400, L(27)80313433200, L(28)80313433200, L(29)2329089562800, L(30)2329089562800, L(31)72201776446800, L(32)144403552893600, L(33)144403552893600, L(34)144403552893600, L(35)144403552893600, L(36)144403552893600, L(37)5342931457063200, L(38)5342931457063200, L(39)5342931457063200, L(40)5342931457063200...可以看到大约在 k42 左右L(k) 就会超过 1e16。实际上经过更精确的计算L(42) 已经远大于 1e16。这意味着在我们最终的求和式S(n) Σ x * (floor(n/L(x-1)) - floor(n/L(x)))中x 只需要从 2 枚举到一个很小的上限比如 60 就绝对安全了因为当 L(x-1) n 时floor(n/L(x-1))为 0后续项对求和没有贡献。因此算法的核心步骤就清晰了预处理出所有需要的 L(k) 值直到 L(k) 1e16。对于给定的 n从 x2 开始枚举计算term x * (floor(n/L(x-1)) - floor(n/L(x)))并累加到答案中。当L(x-1) n时循环终止。这里有一个关键的实现细节L(k) 可能非常大远超 64 位整数范围1e19 量级。在枚举过程中我们计算floor(n / L(x-1))时如果 L(x-1) 本身已经大于 n那么除法的结果就是 0。但如果我们试图计算一个大于 1e19 的 L(k) 值在 C 等语言中可能会发生溢出。因此在预处理 L(k) 序列时我们需要在它即将超过一个很大的界限比如 1e18 或直接与 n 比较时将其设置为一个标记值例如INF n1以表示它已经大于 n这样在后续计算中可以直接判断并跳过。4. 算法实现细节与边界处理理论清晰后我们来看看具体的代码实现。以 C 为例我们需要处理大数运算和溢出问题。首先预处理 L 数组。我们用一个vectorlong long来存储 L(k)但这里long long可能不够最大约9e18。由于 n 最大为 1e16我们可以设定一个上限INF n 1或者一个比 n 大的数如5e18当 L(k) 的计算结果超过这个上限时我们就将其设置为INF表示“无穷大”对于当前的 n 来说。#include bits/stdc.h using namespace std; using ll long long; const ll MOD 1e9 7; const ll INF 5e18; // 一个足够大的数大于 1e16 即可 ll lcm(ll a, ll b) { // 计算 lcm(a, b)需要防止中间结果溢出 // lcm a / gcd(a, b) * b先除后乘 ll g __gcd(a, b); // 判断 (a / g) * b 是否会溢出 if (a / g INF / b) { return INF; // 会溢出返回 INF } return a / g * b; } vectorll precompute_L(int max_k) { vectorll L(max_k 1, 1); L[0] 1; // 方便下标对齐L[1]对应LCM(1) for (int k 1; k max_k; k) { L[k] lcm(L[k-1], k); if (L[k] INF) { L[k] INF; } } return L; }这里max_k需要设多大根据之前的分析L(k) 增长极快k 大约在 40-50 之间就会超过 1e16。为了安全起见我们可以设置max_k 60或100。预处理一次对于所有查询都可以复用。接下来是计算 S(n) 的主函数。我们需要枚举 x注意我们的推导中 x 从 2 开始对应 L(x-1) 和 L(x)。在代码中我们需要访问 L[x-2] 和 L[x-1]如果 L 数组下标从 1 开始存储 LCM(1..k)。ll solve(ll n, const vectorll L) { ll ans 0; // x 从 2 开始枚举对应使用 L[x-2] 和 L[x-1] // L[0] 1 (LCM of empty set or LCM(1)?)我们让 L[1] LCM(1)1, L[2]LCM(1,2)2, ... // 更清晰的做法L[i] 表示 LCM(1,2,...,i) // 则条件为L[x-1] | n 且 L[x] ∤ n // 计数floor(n / L[x-1]) - floor(n / L[x]) for (int x 2; ; x) { // 如果 L[x-1] 已经大于 n则后续所有项都为 0终止循环 if (L[x-1] n) { break; } ll cnt1 n / L[x-1]; ll cnt2 (L[x] n) ? 0 : (n / L[x]); // 如果 L[x] 大于 n则 floor(n/L[x]) 0 ll cnt (cnt1 - cnt2) % MOD; ans (ans (cnt * (x % MOD)) % MOD) % MOD; } return ans; }这里有几个非常重要的细节和边界情况循环终止条件当L[x-1] n时floor(n / L[x-1])为 0并且对于所有更大的 xL[x-1]只会更大因为 L 序列非递减所以后续项都为 0可以安全终止。我们的预处理保证了当 L[k] 超过 INF 时被设为 INF而 INF n所以这个判断是有效的。防止中间溢出在计算cnt1和cnt2时我们直接用n / L[x-1]。由于 n ≤ 1e16L[x-1] 至少为 1这个除法在 64 位整数范围内是安全的。cnt是两个long long的差也在安全范围内。取模操作题目要求结果对 MOD1e97 取模。注意cnt和x都需要先取模再相乘最后累加时也要取模以防止溢出。因为cnt最大可以是 n约1e16x最大约 50乘积约 5e17仍在long long范围内约9e18所以这里不先取模也可以但先取模是更安全的做法。x 的上限理论上循环可能持续到 x 很大但由于 L[x-1] 增长极快实际上循环次数很少大约几十次。所以算法的时间复杂度是 O(K)其中 K 是满足 L(K) ≤ n 的最大索引对于 n1e16K 大约在 40 左右。这是一个常数时间算法效率极高。L[0] 的处理在公式中我们使用了 L(x-1)。当 x2 时需要 L(1)。所以我们的 L 数组最好从索引 1 开始存储 L(1)。我们可以让L[0]无意义或设为 1从L[1]开始是 LCM(1)1。在代码中我倾向于让L[i]代表 LCM(1..i)。那么precompute_L函数需要调整一下循环。让我们修正一下预处理和主逻辑vectorll precompute_L(int max_k) { vectorll L(max_k 2, 1); // 多分配一点索引从1开始到max_k1 L[1] 1; // LCM(1) for (int k 2; k max_k 1; k) { L[k] lcm(L[k-1], k); if (L[k] INF) { L[k] INF; } } return L; } ll solve(ll n, const vectorll L) { ll ans 0; // 枚举 x x 从 2 开始 // 我们需要 L[x-1] 和 L[x] for (int x 2; ; x) { if (L[x-1] n) break; // L[x-1] 对应 LCM(1..x-1) ll cnt1 n / L[x-1]; ll cnt2 (L[x] n) ? 0 : (n / L[x]); ll cnt cnt1 - cnt2; cnt % MOD; ans (ans cnt * (x % MOD)) % MOD; } return ans; }在主函数中我们预处理一个足够大的 L 数组比如precompute_L(60)然后对于每个输入的 n 调用solve(n, L)。5. 容斥原理在此题中的具体体现与思维延伸虽然我们最终的公式S(n) Σ x * (floor(n/L(x-1)) - floor(n/L(x)))看起来简洁但它的推导过程深刻体现了容斥原理的思想。让我们再回顾一下我们想统计有多少个 i 满足f(i) x。这等价于 i 满足两个条件L(x-1) | i属于集合 AL(x) ∤ i不属于集合 B满足条件1的 i 的个数是|A| floor(n / L(x-1))。 这些 i 中不满足条件2的即同时满足条件1和L(x) | i的个数是|A ∩ B| floor(n / L(x))因为L(x)是L(x-1)的倍数所以L(x) | i必然蕴含L(x-1) | i。 因此满足“在A中但不在B中”的 i 的个数就是|A| - |A ∩ B|。这正是容斥原理在两个集合情况下的直接应用|A \ B| |A| - |A ∩ B|。在这道题里容斥原理帮助我们干净利落地从“同时满足两个条件”的正面计数转化为了“满足条件A但不满足条件B”的差值计数从而得到了易于计算的表达式。这种“先计算一个大的集合再减去其中不满足额外条件的子集”的思路在数论计数问题中非常常见。我们可以进一步延伸这个思想。如果问题变得更复杂比如求f(i) ≥ x的 i 的个数或者求f(i)的某种加权和容斥原理依然可以派上用场。本质上我们是在根据L(k)序列定义的“整除层级”对 1 到 n 的数进行分层每一层的数具有相同的f(i)值。容斥原理帮助我们清晰地定义了层与层之间的边界。注意在实现时务必注意数据类型的范围。L(k)在 k 较小时就可能超过 64 位有符号整型的最大值~9.22e18。使用long double进行中间计算比较或者像我们上面那样使用一个自定义的INF阈值来防止溢出是两种常见的策略。在比赛中通常n ≤ 1e16L(50)左右就会超过这个值所以将INF设为1e185或5e18是安全的。如果 n 更大可能需要使用__int128或高精度计算。6. 实战测试与常见问题排查为了确保算法的正确性我们需要用一些测试用例来验证。我们可以写一个暴力程序计算小范围 n比如 n ≤ 10000的 S(n)然后与我们优化后的算法结果进行对比。暴力计算 f(n) 的参考代码ll f_brute(ll n) { for (ll x 2; ; x) { if (n % x ! 0) return x; } } ll S_brute(ll n) { ll sum 0; for (ll i 1; i n; i) { sum f_brute(i); } return sum; }对比测试时要注意取模。我们的优化算法是对 MOD 取模后的结果而暴力计算可能得到真实的和可能非常大。我们需要对比的是S_brute(n) % MOD和solve(n, L)。在测试过程中可能会遇到以下问题结果错误首先检查 L 数组计算是否正确。可以打印出前几十个 L[k] 的值与已知数列1, 2, 6, 12, 60, 60, 420, 840, 2520...进行比对。一个常见的错误是在lcm函数中先乘后除导致溢出即使使用了long long。务必使用a / g * b的顺序并在乘法前判断是否溢出。循环无法终止如果L[x-1]没有正确设置为 INF当它真实值超过long long范围时可能会因为溢出而变成一个很小的数甚至负数导致L[x-1] n的判断永远不成立循环变成死循环。确保你的lcm函数在溢出时返回一个大于 n 的标记值如 INF。答案偏小检查容斥公式是否正确。count(x) floor(n/L(x-1)) - floor(n/L(x))。确保你在累加时用的是x * count(x)并且 x 是从 2 开始。另外检查L数组的下标是否与公式中的 x 对应正确。一个实用的调试方法是对于小的 n手动计算每个 x 对应的count(x)看看是否与暴力统计中f(i)x的 i 的个数一致。取模错误确保在加法和乘法运算的每一步都及时取模防止中间结果溢出。虽然本题中cnt * x可能不会超过 1e18对于 n1e16, x50但养成好习惯很重要。同时注意cnt在减法后可能是负数吗在我们的定义中L(x)是L(x-1)的倍数所以floor(n/L(x-1)) floor(n/L(x))cnt是非负的可以放心。让我们用一个具体的例子来走一遍流程。假设 n 10。预处理 L: L[1]1, L[2]2, L[3]6, L[4]12 (10)后续 L[k] 都 10设为 INF。计算 S(10):x2: cnt floor(10/L[1]) - floor(10/L[2]) 10/1 - 10/2 10-55。贡献: 2*510。x3: cnt floor(10/L[2]) - floor(10/L[3]) 10/2 - 10/6 5-14。贡献: 3*412。x4: cnt floor(10/L[3]) - floor(10/L[4]) 10/6 - 10/12 1-01。贡献: 4*14。x5: L[4]1210循环终止。总和 S 10124 26。暴力验证f(1)2, f(2)3, f(3)2, f(4)3, f(5)2, f(6)4, f(7)2, f(8)3, f(9)2, f(10)3。求和2323242323 26。结果正确。7. 总结与同类问题思路拓展回顾这道“Strange Function”其解题过程是一个典型的将定义转化为数学性质再利用数论工具进行高效计算的范例。关键步骤在于理解与转化认识到 f(n)x 等价于 n 能被 1,2,...,x-1 整除但不能被 x 整除进而转化为与 LCM 序列相关的问题。改变计数视角从逐个计算 f(i) 转为统计有多少个 i 使得 f(i) 等于某个 x。这是优化计数问题的常用技巧。应用容斥原理计算满足“能被 A 整除但不能被 B 整除”的数的个数容斥给出了简洁公式。利用 LCM 序列的快速增长特性这使得需要枚举的 x 非常少将问题复杂度降为常数。这种思路可以拓展到许多类似问题。例如如果函数定义改为“f(n) 是最小的正整数 x使得 x 与 n 互质”那么问题就与欧拉函数或质因子有关。或者如果要求计算g(n) f(1) f(2) ... f(n)其中 f(n) 是 n 的因子个数那就是另一个经典的数论分块问题。在竞赛中遇到这种定义“奇怪”的函数不要被吓到。第一反应应该是尝试枚举小数据观察 f(n) 的输出规律寻找其与 n 的质因数分解、整除性、LCM、GCD 等基本数论概念之间的联系。一旦找到了这种联系问题往往就迎刃而解了。最后在代码实现上处理可能溢出的大数运算需要格外小心。设定安全阈值、使用__int128如果环境支持、或者像本题一样利用问题本身的界限n 有限来避免计算过大的 L(k)都是实用的技巧。把数学推导严谨地翻译成高效且健壮的代码是解决这类题目的最后一步也是至关重要的一步。