1. 项目概述从一道国赛真题看数论与编程的深度结合“矩阵求和”这个题目乍一看像是考察二维数组遍历和基础循环的编程题。但如果你真这么想并且试图用暴力枚举去解决2018年蓝桥杯国赛的这道题那等待你的很可能不是奖牌而是“运行超时”的无情提示。这道题的精妙之处恰恰在于它披着“矩阵求和”的朴实外衣内里却是一道对数论知识和算法优化能力要求极高的综合题。它考察的远不止是写代码更是将数学思维转化为高效算法的能力。简单来说题目会给你一个 n x n 的方阵矩阵中每个位置 (i, j) 的元素值等于 i 和 j 的最大公约数即gcd(i, j)。你的任务是计算这个矩阵中所有元素之和并对一个很大的数通常是 1e97取模。当 n 的值很大时比如达到 1e5 甚至 1e6 的量级直接二重循环计算 gcd 并累加是绝对不可行的。这就需要我们深入挖掘题目背后的数学规律找到那个能将时间复杂度从 O(n²) 降低到 O(n log n) 甚至 O(n) 的“钥匙”。这道题非常适合正在准备算法竞赛如蓝桥杯、ACM的同学尤其是那些已经掌握了基础数据结构希望突破瓶颈学习如何将数论应用于实际编程问题的选手。通过拆解这道题你不仅能学会一道题的解法更能掌握一种“透过现象看本质”将复杂问题转化为已知数学模型并高效求解的思维方式。接下来我们就一层层剥开它的外壳看看核心到底在哪里。2. 核心思路拆解化“矩阵求和”为“数论求和”面对一个 n x n 的矩阵最直观的想法是写两层循环遍历每个位置 (i, j)计算gcd(i, j)并累加。这个算法的时间复杂度是 O(n² * log(min(i, j)))因为每次计算 gcd 也有开销。当 n10⁵ 时运算次数大概是 10¹⁰ 这个量级显然是不可接受的。我们必须寻找更聪明的方法。首先我们明确最终目标是求S(n) Σ(i1 to n) Σ(j1 to n) gcd(i, j)2.1 关键转化引入欧拉函数数论中一个非常经典且强大的工具是欧拉函数 φ(k)它表示小于等于 k 的正整数中与 k 互质的数的个数。这里有一个至关重要的恒等式它是我们解题的基石对于任意正整数 n有 Σ(d|n) φ(d) n。这个公式的意思是n 的所有正因子 d 的欧拉函数 φ(d) 之和等于 n 本身。这个公式的逆用可以让我们用欧拉函数来表示最大公约数。我们考虑 gcd(i, j) 的所有可能取值。设d gcd(i, j)那么 d 必然是 i 和 j 的公约数。反过来对于每一个可能的 d有多少对 (i, j) 满足gcd(i, j) d呢我们可以做如下变换令i d * ij d * j那么条件gcd(i, j) d就等价于gcd(i‘ j’) 1并且1 ≤ i‘ j’ ≤ n/d。于是满足gcd(i, j) d的 (i, j) 对数就等于在1到n/d的范围内满足gcd(i‘ j’) 1的 (i‘ j’) 对数。2.2 推导求和公式现在我们可以重写原始的和式S(n) Σ(i1 to n) Σ(j1 to n) gcd(i, j) Σ(d1 to n) d * [满足 gcd(i, j) d 的 (i, j) 对数]根据上面的变换[满足 gcd(i, j) d 的 (i, j) 对数]等于[满足 gcd(i‘ j’) 1 且 1 ≤ i‘ j’ ≤ n/d 的 (i‘ j’) 对数]。那么在1到m这里m n/d的范围内有多少对 (a, b) 满足gcd(a, b) 1呢这可以通过容斥原理或者再次利用欧拉函数来求。一个更简洁的思路是我们先固定 a看有多少个 b 满足gcd(a, b) 1且1 ≤ b ≤ m。根据定义这个数量就是φ(a)欧拉函数。因此对于所有 a 从 1 到 m满足条件的 b 的数量之和就是Σ(a1 to m) φ(a)。注意这里 (a, b) 是有序对即 (1,2) 和 (2,1) 算作不同的对。所以在1到m的范围内互质有序对的总数就是Σ(a1 to m) Σ(b1 to m) [gcd(a,b)1] Σ(a1 to m) φ(a)不对这里我们重复计算了。实际上对于每一个固定的 a满足条件的 b 的数量是 φ(a)所以总和是Σ(a1 to m) φ(a)。但当我们对 a 求和时已经涵盖了所有有序对。因此[满足 gcd(i‘ j’) 1 且 1 ≤ i‘ j’ ≤ m 的 (i‘ j’) 对数]就等于Σ(a1 to m) φ(a)。注意这里容易产生混淆。更严谨的推导是Σ(a1 to m) Σ(b1 to m) [gcd(a,b)1] Σ(d1 to m) μ(d) * floor(m/d)²这是用莫比乌斯反演的常见形式。而我们通过引入欧拉函数可以得到另一个等价表达式Σ(a1 to m) Σ(b1 to m) [gcd(a,b)1] 2 * Σ(a1 to m) φ(a) - 1。这个公式可以通过计算(a, b)且a ≤ b的互质对数即欧拉函数前缀和与 1 的关系再乘以 2 并减去重复的(1,1)得到。但对于本题的最终求和公式我们可以走一条更直接的经典路径。经典的结论是可以通过莫比乌斯反演或上述有序对思路推导得出Σ(i1 to n) Σ(j1 to n) gcd(i, j) Σ(d1 to n) φ(d) * floor(n/d)²我们来理解一下这个公式对于每一个可能的公约数 d它可能是很多对 (i, j) 的公约数不一定最大。有多少对呢i 和 j 必须都是 d 的倍数。在 1 到 n 的范围内d 的倍数有floor(n/d)个。所以 i 有floor(n/d)种选择j 也有floor(n/d)种选择总共floor(n/d)²对。在这些对中最大公约数可能是 d, 2d, 3d, ...。那么其中最大公约数恰好为 d 的对数是多少呢这正是欧拉函数 φ(d) 可以表示的通过莫比乌斯反演可以证明最大公约数恰好为 d 的 (i, j) 对数为φ(d)。因此所有最大公约数为 d 的数对对总和的贡献是d * φ(d)。但我们的公式是φ(d) * floor(n/d)²这里似乎有点出入。实际上更准确和通用的推导如下 令g(n) Σ(i1 to n) Σ(j1 to n) gcd(i, j)。 我们也可以枚举最大公约数 dg(n) Σ(d1 to n) d * f(n, d)其中f(n, d)是满足gcd(i, j) d且1 ≤ i, j ≤ n的 (i, j) 对数。 令i d*a,j d*b则条件变为1 ≤ a, b ≤ n/d且gcd(a, b) 1。 所以f(n, d) Σ(a1 to m) Σ(b1 to m) [gcd(a,b)1]其中m floor(n/d)。 而Σ(a1 to m) Σ(b1 to m) [gcd(a,b)1]这个式子有一个非常漂亮的结论Σ(a1 to m) Σ(b1 to m) [gcd(a,b)1] Σ(k1 to m) μ(k) * floor(m/k)²。这是莫比乌斯反演的标准形式。 但还有一个更简洁的、利用欧拉函数前缀和的表达式。通过对称性和欧拉函数的性质可以推导出过程略Σ(a1 to m) Σ(b1 to m) [gcd(a,b)1] 2 * Σ(i1 to m) φ(i) - 1。然而对于本题我们不必直接计算f(n,d)。有一个更高效的变换直接利用狄利克雷卷积的性质。我们回到原始和式g(n) Σ(i1 to n) Σ(j1 to n) gcd(i, j) Σ(i1 to n) Σ(j1 to n) Σ(d|gcd(i,j)) φ(d)。 这一步用到了前面提到的恒等式n Σ(d|n) φ(d)将gcd(i,j)用其所有因子的欧拉函数和表示。 交换求和顺序g(n) Σ(d1 to n) φ(d) * Σ(i1 to n) Σ(j1 to n) [d | i 且 d | j]。[d | i 且 d | j]表示 i 和 j 都是 d 的倍数。在 1 到 n 的范围内d 的倍数 i 有floor(n/d)个同样 j 也有floor(n/d)个。 因此Σ(i1 to n) Σ(j1 to n) [d | i 且 d | j] floor(n/d) * floor(n/d) floor(n/d)²。 于是我们得到了最终的核心公式g(n) Σ(d1 to n) φ(d) * floor(n/d)²这个公式就是本题算法的灵魂。它将一个二维的、与 gcd 相关的求和转化为了一个一维的、与欧拉函数相关的求和。时间复杂度从 O(n²) 降到了 O(n)因为我们需要计算 1 到 n 的欧拉函数并对每个 d 计算一次floor(n/d)。而floor(n/d)的值在 d 较大时是成段不变的这又引出了我们下一步的关键优化数论分块。2.3 思路总结与算法选择所以我们的解题路线图非常清晰目标计算S(n) Σ(i1 to n) Σ(j1 to n) gcd(i, j)。核心转化利用恒等式gcd(i, j) Σ(d|gcd(i,j)) φ(d)和交换求和顺序得到S(n) Σ(d1 to n) φ(d) * floor(n/d)²。任务分解任务一高效求出 1 到 n 所有整数的欧拉函数值 φ(1), φ(2), ..., φ(n)。任务二高效计算求和式Σ(d1 to n) φ(d) * floor(n/d)²。由于floor(n/d)对于连续的 d 取值相同我们可以使用数论分块整除分块来加速将复杂度从 O(n) 降至 O(√n)。基于这个思路我们选择的算法组合是线性筛法求欧拉函数数论分块求和。这是一个非常经典且高效的组合拳。3. 核心模块实现线性筛与数论分块有了理论公式接下来就是如何用代码高效实现。这里有两个技术核心如何快速得到欧拉函数数组以及如何利用数论分块加速求和。3.1 线性筛法求欧拉函数数组求 1 到 n 每个数的欧拉函数如果对每个数单独用公式φ(n) n * Π(1 - 1/p)计算复杂度是 O(n√n)对于 n10⁶ 还能接受但对于更大的 n 或多次查询就不够快了。我们需要O(n)的线性筛法。线性筛欧拉筛不仅可以筛出素数还可以在筛的过程中递推求出欧拉函数。原理如下初始化phi[1] 1。遍历 i 从 2 到 n如果 i 是素数未被标记则phi[i] i - 1因为素数与小于它的所有正整数都互质。遍历已得到的素数列表primes中的每个素数p如果i * p n跳出循环。标记i * p为合数。如果p是i的质因子即i % p 0那么phi[i * p] phi[i] * p。这是因为 i 已经包含了质因子 p根据欧拉函数公式φ(i*p) i*p * Π(1-1/p_i)而 i 和 i*p 的质因子集合相同所以φ(i*p) φ(i) * p。如果p不是i的质因子即i % p ! 0那么phi[i * p] phi[i] * (p - 1)。因为 i 和 p 互质根据积性函数性质φ(i*p) φ(i) * φ(p) φ(i) * (p-1)。如果i % p 0跳出内层循环。这是线性筛的关键保证了每个合数只被其最小质因子标记一次。// C 示例代码线性筛求欧拉函数 #include vector using namespace std; const int MAX_N 1000000; // 根据题目要求调整 int phi[MAX_N 1]; vectorint primes; bool is_composite[MAX_N 1]; void euler_sieve(int n) { phi[1] 1; // 注意通常定义 φ(1)1 for (int i 2; i n; i) { if (!is_composite[i]) { // i是素数 primes.push_back(i); phi[i] i - 1; // 素数的欧拉函数值为 i-1 } for (int p : primes) { if (i * p n) break; is_composite[i * p] true; if (i % p 0) { phi[i * p] phi[i] * p; // p是i的质因子 break; // 保证每个数只被最小质因子筛一次 } else { phi[i * p] phi[i] * (p - 1); // p与i互质 } } } }实操心得线性筛的代码模板需要非常熟练。特别注意phi[1] 1的初始化以及if (i % p 0) break;这一行它是保证 O(n) 复杂度的关键。很多同学在记忆时容易忘记break导致每个合数被多次标记虽然结果可能正确但效率退化成了 O(n log log n) 左右。3.2 数论分块加速求和现在我们有了phi[]数组需要计算S(n) Σ(d1 to n) phi[d] * (n/d)²。直接遍历 d 从 1 到 n 求和复杂度 O(n)。当 n 很大时例如 10¹²虽然本题通常 n 在 10⁶~10⁷ 量级但此技巧很重要仍然很慢。观察floor(n/d)当 d 在某个区间[l, r]内变化时floor(n/d)的值是相同的。我们可以找到这些区间然后一次性处理一个区间内的求和。数论分块原理 对于给定的 n 和左端点 l值k floor(n/l)对应的最大右端点 r 满足floor(n/r) k且floor(n/(r1)) k。这个 r 可以直接计算出来r floor(n / floor(n/l))。因此算法步骤如下初始化ans 0l 1。当l n时计算k n / l整数除法。计算r n / k。此时对于所有d ∈ [l, r]都有floor(n/d) k。我们需要计算区间[l, r]对答案的贡献k² * Σ(dl to r) phi[d]。为了快速得到phi[l] phi[l1] ... phi[r]我们需要在预处理欧拉函数时同时计算出欧拉函数的前缀和数组sum_phi[]其中sum_phi[i] phi[1] ... phi[i]。区间和等于sum_phi[r] - sum_phi[l-1]。贡献为(k * k) % MOD * (sum_phi[r] - sum_phi[l-1]) % MOD。注意取模和防止负数。将贡献加到ans并对 MOD 取模。令l r 1进入下一个区间。循环结束ans即为最终结果。// C 示例代码数论分块求和部分 const long long MOD 1000000007; long long solve(int n) { euler_sieve(n); // 预处理 phi 和 sum_phi // 假设我们已经有了 sum_phi[] 数组sum_phi[i] Σ_{j1}^{i} phi[j] long long ans 0; for (int l 1, r; l n; l r 1) { int k n / l; // floor(n/l) 的值 r n / k; // 与 k 相同的最大右端点 // 计算区间 [l, r] 的贡献 long long segment_sum (sum_phi[r] - sum_phi[l - 1]) % MOD; long long contribution (1LL * k * k) % MOD * segment_sum % MOD; ans (ans contribution) % MOD; } // 处理 ans 可能为负数的情况 ans (ans % MOD MOD) % MOD; return ans; }注意事项这里有几个易错点。第一k * k可能会溢出所以在 C 中要使用1LL * k * k将其提升到long long类型再计算。第二segment_sum通过前缀和相减得到可能为负数所以在最后加 MOD 再取模确保非负。第三取模运算%的优先级与乘除法相同要注意使用括号保证运算顺序或者分步计算。4. 完整代码实现与逐行解析将线性筛和数论分块结合起来我们就可以得到本题的完整高效解法。下面给出一个完整的 C 实现并加上详细注释。#include iostream #include vector using namespace std; const int MAX_N 1000000; // 根据题目数据范围调整通常国赛 n 在 1e6 左右 const long long MOD 1000000007LL; long long phi[MAX_N 5]; // 欧拉函数值数组 long long sum_phi[MAX_N 5]; // 欧拉函数前缀和数组 vectorint primes; // 素数表 bool is_composite[MAX_N 5]; // 合数标记 // 线性筛法同时计算欧拉函数 phi[] 及其前缀和 sum_phi[] void init_euler(int n) { phi[1] 1; sum_phi[1] 1; for (int i 2; i n; i) { if (!is_composite[i]) { primes.push_back(i); phi[i] i - 1; // i 是素数 } for (int p : primes) { if (1LL * i * p n) break; // 防止越界 is_composite[i * p] true; if (i % p 0) { phi[i * p] phi[i] * p; // p 是 i 的最小质因子 break; // 关键保证每个数只被最小质因子筛一次 } else { phi[i * p] phi[i] * (p - 1); // p 与 i 互质 } } // 计算前缀和注意取模 sum_phi[i] (sum_phi[i - 1] phi[i]) % MOD; } } // 主求解函数 long long solve(int n) { init_euler(n); // 预处理 long long ans 0; // 数论分块 for (int l 1, r; l n; l r 1) { int k n / l; // floor(n/l) 的值 r n / k; // 当前块的右边界 // 计算区间 [l, r] 的欧拉函数和 long long segment_phi_sum (sum_phi[r] - sum_phi[l - 1]) % MOD; // 计算 (n/l)^2 % MOD long long k_squared (1LL * k * k) % MOD; // 1LL 防止 int 乘法溢出 // 计算当前块的贡献并累加 long long contribution (k_squared * segment_phi_sum) % MOD; ans (ans contribution) % MOD; } // 确保结果非负 return (ans % MOD MOD) % MOD; } int main() { int n; // 假设输入 n 根据题目要求可能有多组测试数据这里以单次为例 // cin n; n 100; // 示例输入 cout solve(n) endl; return 0; }逐行解析与关键点全局定义(MAX_N,MOD)根据题目可能的最大数据范围定义数组大小。MOD 是取模数蓝桥杯这类题常用 1e97。数组与容器phi[]存储欧拉函数sum_phi[]是其前缀和方便区间查询。primes动态存储筛出的素数is_composite[]标记合数。init_euler函数phi[1] 1;是定义。外层循环for (int i 2; i n; i)遍历每个数。if (!is_composite[i])判断 i 为素数初始化phi[i] i - 1并加入素数表。内层循环for (int p : primes)用已得素数筛去合数。if (1LL * i * p n) break;是防溢出和越界的标准写法。if (i % p 0) ... break;是线性筛的核心逻辑决定了欧拉函数的递推公式和筛法的线性复杂度。循环结束后计算sum_phi[i]注意取模。solve函数首先调用init_euler(n)进行预处理。for (int l 1, r; l n; l r 1)是数论分块的经典循环结构。int k n / l;和r n / k;确定了当前块的范围[l, r]其中floor(n/d)都等于k。segment_phi_sum利用前缀和数组O(1)计算出区间[l, r]的φ(d)之和。k_squared计算k² mod MOD使用1LL转换防止int乘法溢出。contribution是当前块对总和的贡献。循环结束后对ans进行最终取模处理确保返回非负数。复杂度分析init_euler是O(n)。solve中的数论分块循环次数约为2√n级别因为l的取值是1, r₁1, r₂1, ...增长很快。整体复杂度为O(n √n)对于 n ≤ 10⁷ 可以轻松应对。5. 常见问题、调试技巧与扩展思考即使理解了算法在实现时也可能遇到各种问题。这里总结一些常见的“坑”和调试技巧。5.1 常见问题与解决结果错误可能是负数原因在取模运算中(a - b) % MOD当a b时结果为负数。我们在计算segment_phi_sum (sum_phi[r] - sum_phi[l-1]) % MOD和最后返回ans时都可能遇到。解决在相减后加上 MOD 再取模。标准写法是((a - b) % MOD MOD) % MOD。上面的代码在最后统一处理了ans但更安全的做法是在计算segment_phi_sum时就处理(sum_phi[r] - sum_phi[l-1] MOD) % MOD。运行超时但 n 并不大原因最可能的原因是线性筛写错了比如漏掉了if (i % p 0) break;这一行。这会导致每个合数被其所有质因子重复标记虽然结果正确但内层循环次数大增复杂度退化。检查仔细核对线性筛的模板特别是break条件。可以写一个小程序输出前20个 phi 值与手算结果对比φ(1)1 φ(2)1 φ(3)2 φ(4)2 φ(5)4 φ(6)2 ...。答案对不上差一点检查 φ(1)确保phi[1]初始化为 1。有些资料或实现中 φ(1) 定义为 0但在此公式推导中必须为 1。检查求和范围公式是Σ(d1 to n)循环时l从 1 开始。确认没有漏掉 d1 的情况。检查取模确认所有乘法和加法操作都及时取模防止中间结果溢出。尤其是k * k在 C 中两个int相乘可能溢出必须用1LL * k * k。验证小数据用暴力算法双重循环计算 n 较小比如 n10, 20时的结果与你的优化算法结果对比这是最有效的调试方法。内存超限原因MAX_N设置得过大或者使用了vectorbool以外的布尔数组每个bool占 1 字节。对于 n10⁷phi、sum_philong long数组约占 160MBis_compositebool约占 10MB在蓝桥杯环境可能接近极限。优化使用bitset或vectorbool来标记合数可以节省 8 倍内存vectorbool是特化的每个元素占 1 bit。如果 n 真的非常大如 10⁸线性筛的内存可能无法承受需要考虑使用时间复杂度稍高但内存更小的筛法如分段筛或者直接使用O(n log n)的普通筛法求欧拉函数。5.2 调试技巧与测试用例单元测试将init_euler和solve函数分开测试。测试init_euler输出前 20 个 phi 值与已知值对比。测试solve编写一个暴力函数long long brute_force(int n)用双重循环计算小 n 的答案。用assert(solve(n) brute_force(n))进行验证。边界测试n 1矩阵只有元素 (1,1)gcd(1,1)1答案应为 1。n 2矩阵为 [[1,1], [1,2]]和为 11125。用公式算d1时 φ(1)floor(2/1)²144d2时 φ(2)floor(2/2)²111总和为5。n 100 或 1000与暴力结果对比。性能测试用n 1e6测试运行时间。在 OJ 环境下O(n) 的算法应该在 1 秒内完成。5.3 扩展思考与变式这道题的本质是计算Σ gcd(i,j)。掌握了这个方法可以解决一系列变式问题二维前缀和形式求矩阵中某个子矩阵[x1,y1]到[x2,y2]的 gcd 和。可以利用容斥原理转化为四个S(n)形式问题的组合但需要注意下标从 1 开始和公式的适配。最大公约数不为 1 的对数求1 ≤ i, j ≤ n且gcd(i,j) 1的 (i,j) 对数。可以用总对数 n² 减去互质对数Σ φ(d) * floor(n/d)²其中 d1 的项就是互质对数需要仔细推导。更直接的是用莫比乌斯反演。最小公倍数求和求Σ(i1 to n) Σ(j1 to n) lcm(i, j)。这比 gcd 求和更复杂需要利用lcm(i,j) i*j / gcd(i,j)然后推导出与 gcd 相关的式子通常也会用到欧拉函数或莫比乌斯反演。多维扩展求Σ gcd(i,j,k)三维。思路类似公式会变为Σ φ(d) * floor(n/d)³。个人体会这道“矩阵求和”题是连接基础编程和中级数论算法的绝佳桥梁。它告诉我们在算法竞赛中看到数据范围大的题目暴力枚举往往是死路。必须静下心来分析问题的数学本质。欧拉函数、前缀和、数论分块整除分块这些工具单独看都不难但组合起来就能解决看似复杂的问题。多积累这样的“组合技”并理解其背后的推导过程比死记硬背模板要重要得多。在平时练习时不妨多问自己这个公式是怎么来的为什么这样优化就快了还有没有其他方法只有这样才能在下一次遇到披着“矩阵”外衣的“数论”题时一眼看穿它的真面目。