莫比乌斯反演与Min_25筛：解决算法竞赛中的高难度数论求和问题-北京尧图网络

1. 项目概述一场算法竞赛中的“硬骨头”如果你参加过算法竞赛尤其是像全国大学生程序设计竞赛ICPC/CCPC或蓝桥杯这类对算法深度和代码效率要求极高的比赛那你一定对“卡常”这个词不陌生。它意味着你的算法思路完全正确复杂度分析也达标但就是差那么一点点无法在规定的时间和内存限制内通过所有测试点。这时候就需要你化身“代码调优师”在算法框架内进行极致的优化。今天要聊的这个项目就是这样一个典型的案例一个结合了莫比乌斯反演和Min_25筛的高难度数论求和问题并且最终需要“卡常”才能通过。这个问题的核心是计算一个与数论函数相关的和式。莫比乌斯反演是处理这类和式的经典工具它能将复杂的、包含整除条件的求和转化为相对简单的形式。然而转化后的新和式往往需要对一个数论函数进行快速的前缀和计算当数据范围巨大比如 n 在 1e10 这个量级时传统的线性筛法就完全无能为力了。这时Min_25筛就闪亮登场了。它是一种亚线性筛法专门用来在极短时间内计算积性函数的前缀和。所以这个项目的技术栈非常清晰理论推导靠莫比乌斯反演高效计算靠Min_25筛最后通过“卡常”技巧将性能压榨到极限。这几乎代表了算法竞赛中数论问题的最高难度梯队——不仅要求参赛者有扎实的数学功底能完成复杂的公式推导还要求有高超的工程实现能力能将理论算法用代码高效、稳定地实现出来并针对评测环境进行微调。接下来我们就一步步拆解这块“硬骨头”。2. 核心思路与数学模型构建面对一个复杂的求和问题直接暴力计算是行不通的。我们的第一步永远是进行数学上的化简与转化寻找可计算的规律。2.1 问题原貌与莫比乌斯反演的应用通常这类问题的原始形式可能类似于求∑_{i1}^{n} ∑_{j1}^{m} f(gcd(i, j))或者更复杂的嵌套形式。其中f是一个定义在正整数上的函数。gcd最大公约数的存在使得求和项之间相互耦合直接计算复杂度是 O(nm)不可接受。莫比乌斯反演的核心思想就在这里发挥作用。它利用莫比乌斯函数μ(n)的性质提供了一个将“整除求和”与“等于求和”相互转换的公式。最常用的一个形式是[n 1] ∑_{d|n} μ(d)其中[ ]是艾弗森括号当括号内条件为真时值为1否则为0。通过设置n gcd(i, j)我们可以将条件gcd(i, j) k转化为[gcd(i, j) k] [gcd(i/k, j/k) 1]进而应用上述公式。经过一系列标准的推导步骤设i i/k,j j/k原问题通常可以转化为如下形式Ans ∑_{k1} f(k) * S(floor(n/k)) * S(floor(m/k))其中S(x) ∑_{d1}^{x} μ(d) * floor(x/d)^2这里以二维求和为例具体形式随问题变化。这个转化是巨大的进步。我们将一个二维的、耦合的求和变成了一个一维的、关于k的求和。并且求和项中包含floor(n/k)和floor(m/k)这意味着当k变大时floor(n/k)的取值会保持不变形成许多“整块”。这提示我们可以使用数论分块也叫除法分块来加速计算将复杂度从 O(n) 降低到 O(√n)。2.2 新瓶颈积性函数前缀和的快速计算数论分块要求我们能快速计算S(x)在任意x处的值。S(x)是另一个和式其核心是计算莫比乌斯函数μ(n)的前缀和M(x) ∑_{i1}^{x} μ(i)。对于x在 1e7 以内我们可以用线性筛预处理出所有μ(i)然后求前缀和查询就是 O(1) 的。但是竞赛题目的数据范围往往故意卡在这个界限之上比如n最大为 1e10。此时我们无法筛出 1e10 以内的所有μ(i)无论是时间还是空间都不允许。这就是引入Min_25筛的原因。Min_25筛最初被设计用来计算积性函数f(n)的前缀和∑_{i1}^{n} f(i)它只需要 O(n^{3/4} / log n) 的时间复杂度和 O(√n) 的空间复杂度。对于n1e10这个复杂度是完全可接受的。虽然μ(n)本身不是积性函数吗它是的。莫比乌斯函数μ(n)是一个经典的积性函数。因此我们可以直接使用 Min_25筛来计算其前缀和M(x)。至此我们的整体算法框架就确定了推导使用莫比乌斯反演将原问题转化为需要计算M(x)的问题。分块对转化后的一维求和式使用数论分块。计算在数论分块需要求M(x)时调用 Min_25筛算法进行计算。优化对 Min_25筛的实现进行极致优化卡常以满足时限要求。注意这里为了叙述清晰以计算M(x)为例。实际题目中S(x)可能更复杂可能需要计算多个积性函数组合的前缀和但核心思想不变——用 Min_25筛解决积性函数在大范围下的前缀和查询问题。3. Min_25筛算法原理深度解析Min_25筛是本题的核心引擎理解其原理对于实现和优化至关重要。它巧妙地将前缀和计算分为两个步骤。3.1 第一步求“质数部分”的和定义g(n, j)表示在 1 到 n 的所有整数中最小质因子大于第 j 个质数P_j的那些数的函数值之和。这里我们用一个“完全积性函数”f(i)来近似我们目标函数f(i)。对于计算μ(i)的前缀和一个常见的选择是令f(i) 1因为对于质数 pμ(p) -1但第一步我们通常先计算一个辅助函数。更常见的做法是我们构造两个函数g0(n, j)表示范围内所有满足条件的数 i 的[i是质数]的贡献计数这里我们用函数f0(i)1来近似。g1(n, j)表示范围内所有满足条件的数 i 的[i是质数] * i的贡献这里我们用函数f1(i)i来近似。初始时g0(n, 0) n - 1减去1因为1不是质数g1(n, 0) n(n1)/2 - 1。然后我们进行状态转移从j-1推到jg(n, j) g(n, j-1) - f(P_j) * [ g(floor(n/P_j), j-1) - g(P_j-1, j-1) ]这个公式的意义是从g(n, j-1)最小质因子大于P_{j-1}的数中减去那些最小质因子恰好等于P_j的数。这些数可以表示为P_j * t其中t的最小质因子大于等于P_j这样才能保证P_j是最小的。所以t的范围是floor(n/P_j)且t满足g(floor(n/P_j), j-1)的条件。但g(P_j-1, j-1)代表的是所有小于P_j的质数构成的集合这部分t是质数但小于P_j其最小质因子就是它本身小于P_j不满足“最小质因子大于等于P_j”的条件所以要加回来。这个过程一直进行到P_j * P_j n为止。最终g0(n, j_max)就近似等于 n 以内的质数个数g1(n, j_max)近似等于 n 以内质数的和。但注意这里我们只是用完全积性函数f得到了一个“近似”它包含了所有质数和一部分合数那些所有质因子都大于P_{j_max}的合数即“大质数”合成的合数。第二步会处理这些合数。3.2 第二步求完整前缀和定义S(n, j)表示在 1 到 n 的所有整数中最小质因子大于等于第 j 个质数P_j的那些数的真实目标函数f(i)值之和。我们最终要求的就是S(n, 1) f(1)。f(1)通常根据函数定义单独处理。S(n, j)可以通过递归计算S(n, j) g(n, j_max) - sp[j-1] ∑_{kj}^{P_k^2 n} ∑_{e1}^{P_k^{e1} n} f(P_k^e) * S(floor(n / P_k^e), k1) f(P_k^{e1})我们来拆解这个公式g(n, j_max) - sp[j-1]这是“质数部分”的贡献。g(n, j_max)是第一步得到的所有“最小质因子大于P_{j_max}”的数的近似和用f算的但我们需要的是真实质数的和sp[...]。sp[j-1]是前j-1个质数的真实函数值之和。所以这部分计算了所有大于等于P_j的质数的贡献。求和符号部分这是“合数部分”的贡献。我们枚举最小质因子P_k从第 j 个开始再枚举这个质因子的指数e。那么一个形如P_k^e * t的数其函数值f(P_k^e * t) f(P_k^e) * f(t)因为f是积性函数且P_k^e与t互质因为t的最小质因子大于P_k。所以贡献是f(P_k^e) * S(floor(n / P_k^e), k1)。最后的 f(P_k^{e1})是补上形如P_k^{e1}的纯质数幂的贡献因为它对应的t1在S(...)中可能未被计入S定义要求最小质因子大于等于P_{k1}对于t1不成立。递归的边界是n P_j时S(n, j) 0。3.3 算法实现中的关键技巧离散化与存储第一步中g(n, j)的n参数会取到所有floor(n / i)i从1到n的值这些值只有 O(√n) 个不同的。我们可以用两个数组id1[ ]和id2[ ]来分别存储x √n和x √n的离散化下标从而将空间复杂度降到 O(√n)。递归与记忆化第二步的S(n, j)是递归计算的且会有大量重复状态。必须用哈希表或数组进行记忆化搜索否则复杂度会退化。预处理需要预处理出一定范围内如 √n 以内的质数列表P[]、质数的真实函数前缀和sp[]以及通过线性筛得到的f(i)在小范围内的值用于递归边界和小范围直接计算。实操心得Min_25筛的代码实现有固定的“板子”但理解和记忆每个数组的含义至关重要。g0/g1、id1/id2、sp0/sp1这些变量在代码中频繁出现建议在写代码时用详细的注释标明每个变量的含义否则调试起来会非常痛苦。尤其是离散化下标映射很容易搞错。4. 从理论到实践完整实现与卡常优化有了理论武器我们需要把它变成高效的C代码。这里以计算∑_{i1}^{n} μ(i)为例展示一个相对完整的Min_25筛实现框架并融入关键的“卡常”技巧。4.1 基础实现框架首先定义常量和全局变量。假设n最大为 1e10。#include bits/stdc.h using namespace std; using ll long long; using ull unsigned long long; const int MAX_SQRT_N 100000; // 因为 n1e10, sqrt(n)1e5 ll n; int sqrt_n, prime_cnt; int primes[MAX_SQRT_N]; bool is_prime[MAX_SQRT_N]; ll sp0[MAX_SQRT_N]; // sp0[i] sum_{j1}^{i} f0(prime[j]) 对于μ(n)质数处f0(p) -1 ll sp1[MAX_SQRT_N]; // 可能不需要这里为了通用性保留 // 离散化相关 ll w[MAX_SQRT_N * 2]; // 存储所有不同的 n/i 值 int w_cnt; int id1[MAX_SQRT_N], id2[MAX_SQRT_N]; // 映射下标 ll g0[MAX_SQRT_N * 2]; // 第一步的 g0 数组 // 记忆化 S(n, j) 的结果 unordered_mapll, ll memo_S;第一步线性筛预处理void linear_sieve(int limit) { fill(is_prime, is_prime limit 1, true); is_prime[0] is_prime[1] false; prime_cnt 0; sp0[0] 0; for (int i 2; i limit; i) { if (is_prime[i]) { primes[prime_cnt] i; // 计算 sp对于 μ(p) -1 sp0[prime_cnt] sp0[prime_cnt - 1] (-1); // sp0 记录质数处的 f0 和 } for (int j 1; j prime_cnt i * primes[j] limit; j) { is_prime[i * primes[j]] false; if (i % primes[j] 0) break; } } }第一步Min_25筛第一部分构造g数组void min25_step1() { w_cnt 0; // 离散化所有可能的 n/i for (ll l 1, r; l n; l r 1) { r n / (n / l); w[w_cnt] n / l; // 初始化 g0 f0(i)1所以 g0初始值是 (n/l) - 1 去掉1 if (w[w_cnt] sqrt_n) id1[w[w_cnt]] w_cnt; else id2[n / w[w_cnt]] w_cnt; g0[w_cnt] w[w_cnt] - 1; // 对应公式中的 g(n, 0) } // DP 更新 g0 for (int j 1; j prime_cnt; j) { ll p (ll)primes[j]; if (p * p n) break; for (int i 1; i w_cnt w[i] p * p; i) { ll x w[i] / p; int idx (x sqrt_n) ? id1[x] : id2[n / x]; // 状态转移: g(n, j) g(n, j-1) - f(p) * [g(n/p, j-1) - g(p-1, j-1)] // 对于 f0(i)1, f(p)1。 g(p-1, j-1) sp0[j-1] g0[i] - (g0[idx] - sp0[j - 1]); } } }第二步递归计算 S(n, j)ll S(ll x, int j) { if (x primes[j] || x 1) return 0; // 记忆化 int idx (x sqrt_n) ? id1[x] : id2[n / x]; ll key (ll)idx * prime_cnt j; // 构造一个唯一的key if (memo_S.count(key)) return memo_S[key]; // 质数部分贡献 g(x, j_max) - sp[j-1] // g0[idx] 对应的是用 f0 近似的“质数”计数我们需要的是真实质数的 μ 和即 - (质数个数) // 所以这部分贡献是 ( - g0[idx] ) - ( - sp0[j-1] ) sp0[j-1] - g0[idx] ll ret sp0[j - 1] - g0[idx]; // 合数部分贡献 for (int k j; k prime_cnt; k) { ll p primes[k]; if (p * p x) break; ll pe p; // 枚举指数 e for (int e 1; pe * p x; e) { ll pe_next pe * p; // μ(p^e) 的值 e1时为-1 e2时为0 ll f_pe (e 1) ? -1 : 0; if (f_pe ! 0) { ret f_pe * S(x / pe, k 1); } // 加上 f(p^{e1}) ll f_pe_next (e 1 1) ? -1 : 0; ret f_pe_next; pe pe_next; } } memo_S[key] ret; return ret; } // 主函数计算 M(n) ∑_{i1}^{n} μ(i) ll calc_mu_prefix(ll x) { n x; sqrt_n sqrt(n); linear_sieve(sqrt_n); min25_step1(); memo_S.clear(); // S(n, 1) 计算的是最小质因子P1的数的和不包括1。 μ(1)1。 return S(n, 1) 1; }4.2 “卡常”优化技巧实录上面的代码框架是正确的但在n1e10时很可能超时。我们需要进行一系列优化。递归剪枝在S(n, j)函数中合数部分枚举k时条件p * p x就break这已经是一个重要剪枝。此外当f(P_k^e) 0时对于μ函数e2时即为0内层循环可以直接break因为更高次幂的贡献也是0。记忆化优化使用unordered_map可能成为性能瓶颈。因为n/i的值只有 O(√n) 个j的范围是质数个数 O(√n / log n)。我们可以用一个二维数组mem[SQRT_N][PRIME_CNT]来存储但PRIME_CNT可能很大~1e4二维数组太大。更常用的技巧是利用x的离散化下标idx和j构造一个一维键值。// 将 key idx * prime_cnt j 改为 // 由于 j 较小我们可以用 idx * (prime_cnt1) j或者直接用 idx 20 | j 如果质数个数小于2^20 // 但更高效的是使用自定义哈希表或直接开一个大数组用 -1 初始化。 vectorll memo_table(w_cnt * (prime_cnt 2), LLONG_MIN); ll mem memo_table[idx * (prime_cnt 1) j]; if (mem ! LLONG_MIN) return mem; mem ret; return ret;用vector和特殊标记值比unordered_map快很多。整数除法与循环优化将ll替换为unsigned long long在某些情况下能加快除法速度。在min25_step1的双重循环中内层循环的w[i] p * p判断可以提前break。避免在循环内进行不必要的类型转换。预处理与缓存sp0[j-1]在循环中频繁使用可以提前存到局部变量。primes[k]也可以提前取出。递归函数S的内联与参数传递将S函数声明为inline并且将频繁使用的n作为全局变量避免参数传递。将primes、sp0等数组也作为全局变量访问。编译器优化使用#pragma GCC optimize(O3)和#pragma GCC target(avx2)等指令注意竞赛环境是否允许。使用-Ofast编译选项本地测试。将递归函数改为迭代对于 Min_25 第二步递归结构清晰改迭代复杂通常优先优化递归本身。针对μ(n)的特殊优化由于μ(n)在质数幂p^e (e2)时为0这使得我们在第二步枚举合数时一旦e2就可以停止当前质数的枚举因为后续更高次幂的贡献都是0。这能显著减少递归调用次数。优化后的S函数核心部分可能如下inline ll S(ll x, int j) { if (x primes[j]) return 0; int idx (x sqrt_n) ? id1[x] : id2[n / x]; int key idx * (prime_cnt 1) j; if (memo[key] ! INF) return memo[key]; // INF 是一个不可能的值如 -1e18 ll ret sp0[j-1] - g0[idx]; // 质数部分贡献 for (int k j; k prime_cnt; k) { ll p primes[k]; ll p2 p * p; if (p2 x) break; // 重要剪枝 ll pe p; // e1 // μ(p) -1 ret (-1) * S(x / pe, k 1); // 对于 μ e1 时已经加了 f(p^{11})μ(p^2)0所以这里不加。 // 并且由于 e2 时 μ(p^e)0直接跳出内层循环 // 因此对于 μ每个质数 k 只枚举 e1 的情况。 // 注意这里需要加上 f(p^{e1}) 当 e1 时即 μ(p^2)0所以实际上没加。 // 但框架保留对于其他函数可能需要。 // break; // 对于μ函数e2无贡献所以直接break外层循环的当前k不对是内层循环只执行一次。 // 更准确地说我们不需要内层 for e 循环了。 // 所以合数部分循环简化为 // ret (-1) * S(x / p, k 1); // 并且不需要再加 f(p^{e1})因为它是0。 } memo[key] ret; return ret; }注意这个简化是针对μ(n)的。对于一般的积性函数仍需保留内层对指数e的循环。踩坑记录最大的一个坑是离散化下标映射。在min25_step1中w[i]存储的是n / l它是递减的。但在S函数中我们需要根据x找到对应的idx。如果x是某个n / l那么它一定在w数组中。映射id1和id2必须正确初始化。一个常见的错误是混淆了id2[n / x]的用法确保n / x是整数且x sqrt_n时n / x sqrt_n这样才能用id2数组索引。务必反复检查这段代码。5. 系统集成与数论分块应用现在我们已经有了计算M(x) ∑ μ(i)的工具。回到最初莫比乌斯反演后的问题。假设我们最终得到了一个需要计算如下形式和的式子Ans ∑_{k1}^{min(n,m)} f(k) * F(floor(n/k)) * G(floor(m/k))其中F(x)和G(x)都可能包含M(x)。数论分块的流程如下ll solve(ll n, ll m) { if (n m) swap(n, m); ll ans 0; for (ll l 1, r; l n; l r 1) { r min(n / (n / l), m / (m / l)); ll n_div n / l; ll m_div m / l; // 假设 F(x) ∑_{i1}^{x} μ(i) * i 我们需要计算 F(n_div) 和 F(m_div) // 这里用我们的 Min_25 筛函数 calc_F_prefix 来计算 ll F_ndiv calc_F_prefix(n_div); ll F_mdiv calc_F_prefix(m_div); // 假设 f(k) 的前缀和可以快速计算记为 sum_f(l, r) ll sum_f_lr sum_f(r) - sum_f(l-1); ans sum_f_lr * F_ndiv * F_mdiv; } return ans; }这里calc_F_prefix需要根据F(x)的具体定义来实现。如果F(x)就是M(x)那直接用calc_mu_prefix即可。如果F(x)是∑ μ(i)*i那么我们需要修改 Min_25筛使其能计算这个函数的前缀和。这需要我们在第一步维护g1近似质数和并在第二步的递归公式中正确处理f(i)μ(i)*i这个函数。性能考量在数论分块中我们会多次调用calc_F_prefix(x)x的取值是所有的floor(n/i)和floor(m/i)总共 O(√n √m) 个不同的值。如果每次调用都重新运行整个 Min_25筛代价是无法承受的。因此我们需要预处理。标准的做法是以最大的n和m作为上界N_max。运行一次 Min_25筛的第一步得到全局的g0,g1等数组。在数论分块循环中对于每个需要查询的x调用第二步的S(x, 1)函数。由于S(x, j)是递归记忆化的所有查询可以共享同一个记忆化表格从而避免重复计算。这意味着我们的代码结构需要稍作调整将 Min_25筛的第一步初始化设为全局一次然后提供一个查询函数query_F(x)内部调用S(x, 1)。6. 常见问题、调试技巧与扩展即使理解了算法实现过程也充满挑战。6.1 常见问题排查清单问题现象可能原因排查方法结果错误小数据对不上暴力解。1. 莫比乌斯反演公式推导错误。2. Min_25筛中f(1)处理不当。3. 离散化下标映射错误。4. 递归边界条件错误。1. 用 n,m 很小如100的数据暴力计算原式和推导式对比。2. 检查S(n,1)返回值是否加了f(1)。3. 输出id1,id2,w数组手动验证几个x的映射是否正确。4. 单步调试S函数查看递归树。结果正确但运行超时。1. 未使用记忆化或记忆化效率低。2. 递归剪枝不充分。3. 常数过大如频繁使用map、%运算。1. 改用数组特殊标记值的记忆化方式。2. 确认对于μ(n)合数部分循环是否利用了μ(p^e)0 (e2)的优化。3. 使用快读减少函数调用使用局部变量。运行时错误段错误。1. 数组开小了。2. 递归过深导致栈溢出。1. 计算w_cnt的最大值它约为2*sqrt(n)。确保相关数组大小足够。2. 评测环境栈空间可能有限。可以考虑将递归改为显式栈较复杂或者尝试优化递归参数减少深度。通常 Min_25 递归深度不会太大。大数据n1e9时结果不对。整数溢出。将所有涉及n * n或大数累加的变量改为long long或unsigned long long。检查p*p是否可能溢出可改为p n / p的判断方式。6.2 调试技巧对拍写一个暴力程序用于小范围数据n,m 1000验证。这是最有效的方法。中间输出在 Min_25 第一步结束后输出g0数组的前若干项与“小于等于 w[i] 的质数个数”的暴力计算结果对比。单元测试单独测试calc_mu_prefix(x)函数与线性筛预处理的结果对比x 在 1e6 以内。静态检查仔细核对所有数组下标确保没有越界。特别是id2[n/x]当x很大时n/x很小要确保这个值在数组范围内。6.3 扩展与总结掌握了“莫比乌斯反演 Min_25筛卡常”这套组合拳你就具备了解决竞赛中绝大多数高难度数论求和问题的能力。这套方法的变体可以用于计算很多积性函数的前缀和例如φ(n)欧拉函数的前缀和。σ_k(n)除数函数的前缀和。λ(n)刘维尔函数的前缀和。以及这些函数的组合形式。其核心思想永远是通过反演化简问题将瓶颈转化为积性函数前缀和问题再用亚线性筛法解决它。最后关于“卡常”我想说它本质上是工程优化是在算法理论复杂度已经最优的前提下针对具体语言、编译器和硬件环境的微调。它不应该成为你思考问题的首要部分。正确的做法是先写出清晰、正确的代码通过小数据测试。只有当算法正确但超时时才像雕刻家一样一点点地打磨代码的性能。优先使用复杂度更优的算法其次才是常数优化。而 Min_25筛本身就是一个通过巧妙设计将常数控制得相当好的算法。理解它掌握它你就能在数论问题的海洋里拥有了一艘强大的破冰船。

莫比乌斯反演与Min_25筛：解决算法竞赛中的高难度数论求和问题

相关推荐

数学建模实战：从数据预处理到模型选择，解析古代玻璃成分分析全流程

AutoPass：基于LLM智能体的编译器性能自动化调优框架

全志 T153 智能家居面板实战：从 USB 触摸屏、MQTT 到 LVGL UI 性能优化

更多资讯

OSM路网入Neo4j：Java路由服务的Dijkstra/A*算路与接口封装

Winform中英文实时切换：XML语言包与控件树递归实践

C++实现大数相乘的算法

CS229机器学习复习主线：从监督学习到SVM与EM的推导链

声源定位精度瓶颈：GCC-PHAT时延估计与MUSIC空间谱实践

YOLO26最新创新改进系列：样条KAN令牌主干，1×1 卷积给每个通道一个固定权重：让 YOLO26 的主干用样条函数去表达非线性！

今日更新

PLFM_RADAR：从被动告警到主动感知的微服务健康监控系统

SPI隔离设计实战：MAX14431BAEE+与STM32F407工业级应用

数字隔离器如何根治SPI通信地环路干扰：MAX14431与TM4C129实战

本周精选

面向复杂光照与角度的宠物面部识别鲁棒性增强策略

基于DRV8818与STM32的双极步进电机驱动方案实战解析

DRV8818PWPR+PIC18F45K80双极步进电机控制方案实战

本月热门

面向复杂光照与角度的宠物面部识别鲁棒性增强策略

基于DRV8818与STM32的双极步进电机驱动方案实战解析

DRV8818PWPR+PIC18F45K80双极步进电机控制方案实战

读完文章，想聊聊您的网站？