1. 项目概述从一道国赛真题看数论与编程的深度结合如果你刷过蓝桥杯的真题尤其是国赛级别的题目一定会对那种融合了算法思维和深厚数学功底的题目印象深刻。今天要拆解的正是2018年蓝桥杯软件类B组国赛的第六题——“矩阵求和”。这道题远不止是一个简单的二维数组累加问题它更像是一个精巧的“伪装者”表面是矩阵计算内核却是一道经典的数论题目核心考察点是对欧拉函数的深刻理解与高效应用。很多选手在赛场上一看题目描述是矩阵可能下意识地就去想二维前缀和或者动态规划结果一头扎进去才发现时间复杂度根本不允许这正是出题人设下的“思维陷阱”。这道题的价值在于它完美地诠释了竞赛编程中一个至关重要的能力将实际问题抽象转化为数学模型并运用合适的数学工具进行高效求解。无论是正在备赛蓝桥杯的选手还是希望提升自己数论与算法结合能力的开发者深入剖析这道题都能带来极大的收获。它不仅教你欧拉函数怎么用更教你什么时候该想到用它。2. 问题本质解析剥开矩阵的“外壳”我们先来还原一下题目的典型描述。给定一个 n x n 的矩阵矩阵中位于第 i 行第 j 列的元素的值定义为gcd(i, j)即整数 i 和 j 的最大公约数。题目要求计算这个矩阵所有元素之和并对最终结果取模通常模数为1e97。用公式表示就是求S(n) Σ_{i1}^{n} Σ_{j1}^{n} gcd(i, j)当 n 较小时比如 n10我们可以直接暴力二重循环计算。但国赛的数据规模 n 通常会非常大可能达到10^6甚至10^7级别。此时 O(n²) 的暴力算法完全不可行时间复杂度爆炸。因此我们必须寻找更优的解法。2.1 核心转化从枚举元素到枚举公约数暴力算法的思路是枚举每一个矩阵元素(i, j)然后计算一次gcd(i, j)。我们需要逆转这个思路。既然矩阵的值是gcd(i, j)那么我们可以考虑矩阵中有多少个元素的值等于某个特定的数 dd 是最大公约数换句话说我们不去计算每对(i, j)的gcd而是去数对于每一个可能的公约数 d有多少对(i, j)满足gcd(i, j) d。记这个数量为cnt(d)。那么矩阵的总和就可以重新表述为S(n) Σ_{d1}^{n} [ d * cnt(d) ]这里d 从 1 枚举到 n。现在问题的关键就变成了如何高效求解cnt(d)。2.2 引入欧拉函数求解cnt(d)的钥匙cnt(d)表示满足1 i, j n且gcd(i, j) d的整数对(i, j)的数量。我们可以做如下变换 令i d * x,j d * y。那么gcd(i, j) d等价于gcd(d*x, d*y) d这进一步等价于gcd(x, y) 1即 x 和 y 互质。同时由于 i 和 j 要小于等于 n所以 x 和 y 必须满足1 x, y n/d。于是问题转化为在1 x, y floor(n/d)的范围内有多少对(x, y)满足gcd(x, y) 1这里的floor(n/d)我们记作m。cnt(d)就等于这个数量。那么如何计算这个数量呢这里就需要欧拉函数 φ(k)出场了。欧拉函数 φ(k) 的定义是小于等于 k 的正整数中与 k 互质的数的个数。一个经典且关键的结论是满足1 x, y m且gcd(x, y) 1的整数对(x, y)的数量等于2 * Σ_{i1}^{m} φ(i) - 1。推导过程如下 我们可以固定 x考虑有多少个 y 满足gcd(x, y) 1。根据定义对于给定的 x满足条件的 y 的数量就是 φ(x)因为 y 在 1 到 m 范围内且与 x 互质。那么对所有 x 从 1 到 m 求和得到Σ_{i1}^{m} φ(i)。但这只是计算了x确定时y的数量即数对了(x, y)中x m的部分。由于(x, y)是无序对gcd(x,y)和gcd(y,x)是同一个元素在矩阵中对应两个位置除非 xy我们实际上需要计算所有有序对(x, y)。 因此总的有序对数量为Σ_{x1}^{m} Σ_{y1}^{m} [gcd(x,y)1]。这个值恰好等于2 * Σ_{i1}^{m} φ(i) - 1。减去 1 是因为当 xy1 时这个互质对在累加中被计算了两次一次作为 (1,1)一次作为对称部分但实际上在矩阵中它只对应一个位置对角线。严谨的证明可以通过莫比乌斯反演得到但上述组合解释更直观。因此我们得到cnt(d) 2 * Σ_{i1}^{m} φ(i) - 1其中m floor(n/d)。注意很多资料和题解会直接给出这个结论。理解其推导过程至关重要这能帮助你在遇到变种题目比如求gcd(i, j)的 k 次方和时知道如何调整思路而不是死记硬背公式。2.3 最终公式与算法框架将cnt(d)代入总和公式S(n) Σ_{d1}^{n} [ d * ( 2 * Σ_{i1}^{floor(n/d)} φ(i) - 1 ) ]至此我们成功将原问题转化为一个需要预处理欧拉函数前缀和然后进行整除分块优化的经典数论问题。算法框架清晰了预处理出1到n的欧拉函数值φ(i)及其前缀和sum_phi[i] Σ_{k1}^{i} φ(k)。遍历d从1到n利用前缀和快速计算sum_phi[floor(n/d)]然后套用公式累加到答案中。但是直接遍历d仍然是 O(n) 的复杂度当 n 很大时如1e7可能依然吃力取决于时限。我们可以利用整除分块或称数论分块进一步优化。3. 核心组件实现欧拉函数与整除分块3.1 欧拉函数的线性筛法计算1到n所有数的欧拉函数最常用的高效算法是线性筛法欧拉筛。它能在 O(n) 的时间复杂度内同时得到素数表和欧拉函数值。原理与代码实现 线性筛的核心是保证每个合数只被其最小的质因子筛掉。在这个过程中我们可以根据质因子关系递推求出 φ(n)。#include vector using namespace std; const int MAXN 1e7 10; // 根据题目数据范围设定 int phi[MAXN]; // 欧拉函数值 vectorint primes; // 保存素数 bool is_prime[MAXN]; // 标记是否为素数可以省略用phi[i]i判断 void euler_sieve(int n) { phi[1] 1; // 定义 for (int i 2; i n; i) { if (!phi[i]) { // 如果phi[i]为0说明i是素数 phi[i] i - 1; // 素数的欧拉函数值为 i-1 primes.push_back(i); } // 遍历已知素数 for (int p : primes) { if (i * p n) break; // 关键此时 p 是 i * p 的最小质因子 if (i % p 0) { phi[i * p] phi[i] * p; // 情况1p是i的质因子 break; // 保证每个数只被最小质因子筛一次 } else { phi[i * p] phi[i] * (p - 1); // 情况2p与i互质 } } } }递推关系解释情况1 (i % p 0): 此时p是i的质因子。那么i * p的质因子与i完全相同。根据欧拉函数公式φ(n) n * Π(1 - 1/p)对于i * p只是n变成了i * p连乘积部分不变。因此φ(i * p) φ(i) * p。情况2 (i % p ! 0): 此时p是i * p的最小质因子且与i互质。根据积性函数性质φ(i * p) φ(i) * φ(p) φ(i) * (p - 1)。实操心得数组phi同时起到了记录函数值和标记访问状态的作用。初始化phi数组为0phi[i]0即表示i未被访问是素数。break语句是线性筛的灵魂它确保了每个合数只被筛一次。当i % p 0时说明p已经是i的最小质因子那么对于后续更大的质数pi * p的最小质因子应该是p而不是p所以应该停止留给后面的i去筛。预处理完成后记得同时计算前缀和数组sum_phi[i] sum_phi[i-1] phi[i]。3.2 整除分块优化求和过程观察最终公式S(n) Σ_{d1}^{n} [ d * ( 2 * sum_phi[floor(n/d)] - 1 ) ]。如果我们直接枚举d复杂度是 O(n)。但注意到当d变化时floor(n/d)的值在很多连续的d区间内是相同的。例如n 10。d从 1 到 10floor(10/d)的值依次为10, 5, 3, 2, 2, 1, 1, 1, 1, 1。 可以看到floor(10/d)2对应了d4,5两个连续的区间floor(10/d)1对应了d6,7,8,9,10一个区间。整除分块就是用来找出这些值相同的连续区间的。对于任意i满足floor(n/i) floor(n/j)的最大j是floor(n / floor(n/i))。这样我们可以将一个区间[i, j]内的贡献一次性计算出来而不是逐个计算。应用在本问题中 设k floor(n / d)。那么对于当前d使得floor(n / d)仍然等于k的d的最大值是j floor(n / k)。区间[d, j]内的所有d其对应的m floor(n/d)都等于k因此公式中的(2 * sum_phi[k] - 1)是常数。 这个区间对总和的贡献就是(2 * sum_phi[k] - 1) * Σ_{xd}^{j} x。 等差数列求和Σ_{xd}^{j} x (d j) * (j - d 1) / 2。优化后的算法步骤预处理phi[1..n]和sum_phi[1..n]。初始化ans 0。令d 1。当d n时循环 a. 计算k n / d。 b. 计算j n / k。 // 整除分块右边界 c. 计算区间和segment_sum (d j) * (j - d 1) / 2。 // 注意取模 d. 计算贡献contribution (2 * sum_phi[k] - 1) * segment_sum。 e. 将contribution累加到ans并取模。 f. 令d j 1进入下一块。输出ans。时间复杂度预处理欧拉筛 O(n)。整除分块部分因为floor(n/d)的值最多只有2√n种所以循环次数为 O(√n)。整体复杂度为 O(n)主要开销在预处理上求和部分极快。4. 完整代码实现与逐行解析下面给出一个完整的 C 实现包含详细的注释并处理了取模运算。#include iostream #include vector using namespace std; typedef long long LL; const int MAXN 10000007; // 根据题目最大数据范围设定例如1e7 const int MOD 1000000007; // 常见的模数 int phi[MAXN]; // 欧拉函数值 LL sum_phi[MAXN]; // 欧拉函数前缀和注意用long long vectorint primes; // 素数表 // 线性筛法预处理欧拉函数及前缀和 void init(int n) { phi[1] 1; for (int i 2; i n; i) { if (!phi[i]) { // i是素数 phi[i] i - 1; primes.push_back(i); } for (int p : primes) { if (i * p n) break; if (i % p 0) { phi[i * p] (LL)phi[i] * p % MOD; // 这里取模不影响phi值本身但防止溢出 break; } else { phi[i * p] (LL)phi[i] * (p - 1) % MOD; } } } // 计算前缀和 sum_phi[0] 0; for (int i 1; i n; i) { sum_phi[i] (sum_phi[i-1] phi[i]) % MOD; } } // 模意义下的整数取逆元用于除法本题未直接使用但常用 LL mod_pow(LL a, LL b) { LL res 1; while (b) { if (b 1) res res * a % MOD; a a * a % MOD; b 1; } return res; } int main() { int n; cin n; // 读入矩阵大小 init(n); // 预处理 LL ans 0; for (int d 1; d n; ) { int k n / d; // floor(n/d) int j n / k; // 整除分块右边界 // 计算 d 到 j 的等差数列和注意取模 // 公式: (首项末项)*项数/2 LL segment_sum ( (LL)(d j) % MOD * (LL)(j - d 1) % MOD ) * mod_pow(2, MOD-2) % MOD; // 使用逆元处理除法 // 更简单且安全的写法在模运算前使用整数运算 // LL segment_sum (LL)(d j) * (j - d 1) / 2; // 由于 (dj) 和 (j-d1) 必有一个是偶数所以整数除法是精确的。 // 但为了通用性下面采用逆元方法实际竞赛中可直接用整数运算。 // 计算当前块的贡献 LL contribution ( (2 * sum_phi[k] % MOD - 1 MOD) % MOD ) * (segment_sum % MOD) % MOD; ans (ans contribution) % MOD; d j 1; // 跳到下一个块 } cout ans endl; return 0; }代码关键点解析数组大小MAXN应比题目可能的最大n稍大防止越界。数据类型sum_phi和中间计算结果使用long long (LL)因为前缀和可能超出int范围。phi数组本身值小于n可以用int。取模运算在计算contribution时2 * sum_phi[k] - 1可能为负数所以先加上MOD再取模确保结果非负(2 * sum_phi[k] % MOD - 1 MOD) % MOD。等差数列求和segment_sum涉及除法/2。在模运算中除法需要转换为乘以2的模逆元。mod_pow(2, MOD-2)利用费马小定理计算了2在模MOD下的逆元。注意因为MOD1e97是质数且2与它互质所以逆元存在。在实际竞赛中如果确定(dj)和(j-d1)中有一个是偶数也可以先进行整数除法再取模这样更快。整除分块循环for (int d 1; d n; )是一个“手动”控制循环变量的写法在循环体内通过d j 1来更新d。5. 常见问题与调试技巧5.1 时间复杂度与空间复杂度权衡问题当n非常大如1e7时预处理phi和sum_phi数组需要约80MB内存int数组约40MBlong long数组约80MB可能接近内存限制。解决方案优化数据类型phi数组元素值不超过n如果n1e7可以用int存储。sum_phi数组元素最大约为n²/2量级当n1e7时约为5e13远超int必须用long long。内存估算1e7个long long是8 * 1e7 ≈ 80MB加上int数组40MB总共120MB。许多竞赛环境内存限制为256MB或512MB这是可以接受的。但如果n更大或内存更紧张需要考虑其他方法。分块处理有一种“杜教筛”的思想可以计算更大n的欧拉函数前缀和但代码复杂。对于蓝桥杯国赛n通常不会大到必须使用杜教筛。5.2 取模运算的陷阱负数的模C中%运算符对负数取模结果是负数。因此在计算(2 * sum_phi[k] - 1) % MOD时如果2*sum_phi[k] - 1是负数结果就是负数。我们需要调整成0到MOD-1之间的数。标准做法是(x % MOD MOD) % MOD。中间结果溢出即使在long long范围内连续的乘法也可能溢出。例如(d j) * (j - d 1)当n很大时j也很大乘积可能超过64位有符号整数范围约9e18。安全做法是在乘法前先取模或者使用__int128如果编译器支持。在竞赛中如果知道n的范围可以评估最大值。对于本题d和j都小于等于n(dj)最大2n(j-d1)最大n乘积最大约2n²。当n1e7时2n²2e14在long long(约9e18) 范围内是安全的。除法与逆元如前所述整除分块中的等差数列求和涉及除以2。最安全的方法是乘以2的模逆元。计算逆元有快速幂法费马小定理要求模数为质数或扩展欧几里得法。5.3 调试与验证技巧小数据暴力对拍编写一个O(n²)的暴力算法用于验证n较小时如n1000优化算法的正确性。这是调试数论题最有效的方法。LL brute_force(int n) { LL res 0; for (int i 1; i n; i) { for (int j 1; j n; j) { res gcd(i, j); // 需要实现gcd函数 } } return res % MOD; }检查欧拉函数值预处理后输出几个小值的phi[i]和sum_phi[i]进行验证。例如phi[1]1, phi[2]1, phi[3]2, phi[4]2, phi[5]4, phi[6]2sum_phi[6] 112242 12检查整除分块区间在循环中打印出每个块的d,j,k值看是否与手动计算一致。使用已知结论验证对于S(n)有一个著名的等式S(n) Σ_{d1}^{n} φ(d) * floor(n/d)²。你可以用这个公式写另一个程序验证结果。这个公式可以通过莫比乌斯反演得到也是本题的另一种推导思路。5.4 性能优化点素数筛的优化线性筛中的内层循环for (int p : primes)可以改为用数组下标访问避免vector的迭代器开销但差别不大。整除分块循环中的计算等差数列求和(dj)*(j-d1)/2可以确保整除因为(dj)和(j-d1)必为一奇一偶。因此可以直接使用整数运算避免模逆元的计算开销这是重要的优化。LL segment_sum (LL)(d j) * (j - d 1) / 2 % MOD;输入输出优化当n很大时使用cin/cout可能较慢。可以加入ios::sync_with_stdio(false); cin.tie(0);来关闭同步流或使用scanf/printf。6. 举一反三变种与拓展思考这道题是数论在竞赛编程中的一个典型应用。掌握其核心思想后可以解决一系列类似问题。6.1 变种一求Σ_{i1}^{n} Σ_{j1}^{n} gcd(i, j)^k即最大公约数的 k 次方和。思路完全一样只是最终公式变为S_k(n) Σ_{d1}^{n} [ d^k * cnt(d) ] Σ_{d1}^{n} [ d^k * (2 * sum_phi[floor(n/d)] - 1) ]同样可以用整除分块优化。需要预处理d^k的前缀和或者在线计算d^k的区间和可能需要快速幂。6.2 变种二求Σ_{i1}^{n} Σ_{j1}^{m} gcd(i, j)(n ≠ m)矩阵不是方阵的情况。推导过程类似S(n, m) Σ_{d1}^{min(n,m)} [ d * cnt(d) ]其中cnt(d)是满足1in, 1jm且gcd(i,j)d的对数。通过变换id*x, jd*y得到cnt(d)等于满足1xn/d, 1ym/d且gcd(x,y)1的对数。这个数量没有像方阵那样简洁的公式但可以通过莫比乌斯反演得到cnt(d) Σ_{k1}^{min(n/d, m/d)} μ(k) * floor(n/(d*k)) * floor(m/(d*k))其中μ(k)是莫比乌斯函数。然后通过两次整除分块或称为数论分块套分块来优化求和。难度显著增加。6.3 拓展思考与其他知识点的联系莫比乌斯反演本题的推导过程实际上是莫比乌斯反演的一个特例。cnt(d)可以通过反演公式cnt(d) Σ_{k1}^{floor(n/d)} μ(k) * floor(n/(d*k))²得到结合Σ_{d|n} φ(d) n这个性质可以推导出我们使用的公式。学习莫比乌斯反演能让你从更高视角理解这类问题。积性函数与线性筛欧拉函数φ(n)是积性函数。线性筛法不仅可以筛素数、求欧拉函数还可以求莫比乌斯函数μ(n)、约数个数d(n)、约数和σ(n)等积性函数模板非常相似。掌握这个模板是数论编程的基础。整除分块的本质它优化的是形如Σ_{i1}^{n} floor(n/i) * f(i)的求和其中f(i)的前缀和可以快速计算。关键在于floor(n/i)的值成块状分布。7. 竞赛实战策略与心得在蓝桥杯等竞赛中遇到此类题目我的策略通常是冷静分析识别模型看到“矩阵”、“gcd求和”、数据规模大要立刻联想到数论转化而不是蛮力模拟。先尝试小规模暴力寻找规律验证猜想。推导公式写在草稿纸上将S(n)转化为Σ d * cnt(d)然后集中精力攻克cnt(d)。如果一时想不起欧拉函数结论可以尝试用莫比乌斯函数推导或者通过小数据打表观察cnt(d)的规律。选择可靠模板线性筛欧拉函数、整除分块都是固定套路赛前必须做到熟练默写避免现场调试。注意取模、负数处理、整数溢出等细节。测试与对拍即使推导无误实现时也可能出错。务必编写暴力程序对拍小数据n100。测试几个典型值如 n1, 2, 10, 100。时间与空间评估根据数据范围n估算预处理数组大小和内存占用。如果n达到1e7使用long long前缀和数组是必要的。确保整体复杂度在 O(n) 以内。个人踩坑记录我曾因为忘记处理公式中2 * sum_phi[k] - 1可能为负数的情况导致结果错误调试了很久。还有一个易错点是在整除分块计算等差数列和时直接写(dj)*(j-d1)/2 % MOD在取模前做除法由于C中整数除法的特性虽然数学上可整除但取模运算介入的顺序错误也会导致问题。后来我统一先对乘法部分取模再用逆元处理除法或者确认在取模前完成所有整数运算才保证了正确性。这道“矩阵求和”题就像一把钥匙打开了利用数论知识优化算法复杂度的大门。理解它不仅能解决一道题更能获得解决一类题的能力。在编程竞赛中这种从具体问题抽象出数学模型并运用高效算法实现的能力才是区分普通选手和高手的关键所在。