P7112【模板】行列式求值：任意模数下的辗转相除消元法

P7112 这道题在 OI/ACM 圈子里的地位有点像高精度、FFT、“平衡树模板”那种人手一份是必须的。题目全名叫“【模板】行列式求值”字面意思就是让你求一个 n×n 方阵的行列式对 p 取模的结果。听起来很线性代数其实竞赛里的行列式和你在线代课上背的定义差别不小关键在于两条数据范围很大n 经常到 600以及 p 不保证是质数。第一点决定了你必须写高效消元而不是按定义展开第二点让很多“背过高斯消元 逆元”的选手直接翻车。这篇文章就展开讲清楚为什么普通高斯消元在这里不够用、辗转相除消元法的原理是什么、最终可复现的模板长什么样以及我在各种题目里复用它时踩过的坑。无论你是刚接触线性代数竞赛题的初学者还是想找个能直接背下来的模板的选手这篇都适用。1. 这道模板题到底在考什么1.1 模板题对选手的真实要求模板题在 OI 圈里指的是题目没有花哨包装考的就是某个算法或数据结构的标准实现给你在考场上“默写”用的。P7112【模板】行列式求值就是这样一张硬核裸题。你拿到手就是一个矩阵和一个模数 p任务只有一个输出行列式模 p。但它的阴险之处在于它强制要求你处理“一般模数”。模数 p 不保证是质数甚至可能是个很小的合数比如 6、8、9。这意味着很多选手熟记的“高斯消元 费马小定理求逆元”解法在 p 是质数时能跑得飞快在这道题里直接失效。所以说它真正考核的并不是“你会不会写高斯消元”而是“能不能对付模数不是质数的情况”。这类模板题的价值在于把算法最核心的部分剥出来逼迫你掌握原理而不是背板子。我见过不少同学可以把质数版本的行列式模板默写得很溜但换到 P7112 就开始乱改最常见的改法是把逆元运算改成扩展欧几里得求逆元或者干脆用浮点数。前者合数照样求不出逆元后者在竞赛环境下必爆精度。真正稳的做法是下面要讲的“不除法”消元。1.2 为什么刻意让模数“不保证是质数”有人会问现实题目里基本都给质数模数为啥这道题要故意整一个非质数因为在很多高级应用里模数确实是任意的。最典型的是矩阵树定理求生成树数量当计数对象本身是一个整数题目如果要求“对某个合数取模”你不能回避。还有生成函数、DP 里的行列式加速、甚至某些数论问题都可能出现模数是合数的情况。竞赛命题人设计 P7112就是在传递一个信号行列式求值模板不应该只建立在“逆元存在”这个狭隘前提下而应该回到线性代数最基本的初等变换性质上去。这才是“模板题”真正该有的样子。下面这张表可以直观看出区别模数情况高斯消元 逆元辗转相除消元法p 为质数可行需要 powmod 求逆元可行不需要逆元p 为合数主元与 p 不互质时失效恒可行无逆元依赖实现难度较低套路清晰中等需要理解交换和符号2. 算法原理从高斯消元到辗转相除2.1 高斯消元求行列式的三板斧行列式计算在竞赛里不靠定义展开那是 O(n!) 的灾难。我们真正依赖的是行列式的三条初等变换规则交换两行行列式变号某一行乘以常数 k行列式乘以 k某一行加上另一行的 k 倍行列式不变。高斯消元用这三条规则把矩阵变成上三角形式然后直接把对角线元素乘起来。因为三角形矩阵的行列式就是主对角线乘积这个结论在线性代数里非常基础。具体来说为了把主元 a[i][i] 下方的 a[j][i] 清成 0我们做行变换a[j] a[j] - (a[j][i] / a[i][i]) * a[i]这里的核心动作就是“除以主元 a[i][i]”。如果你在实数域里做除法毫无压力在模质数 p 的域里做只要主元不是 0就能用费马小定理把倒数求出来所以也毫无压力。这也是为什么大多数教材和模板默认质数模数代码里动不动就是 inv powmod(a[i][i], p-2)。2.2 非质数模数让“除以主元”失效现在假设 p 6主元 a[i][i] 2。请问 2 的逆元存在吗不存在。因为 gcd(2, 6) 2 ≠ 1模 6 的乘法群里压根没有 2 的倒数。于是你没办法写“除以主元”高斯消元的第一板斧就失灵了。但行列式本身是一个确定的整数。比如矩阵 [[1, 2], [3, 4]] 的行列式是 1×4 - 2×3 -2按模 6 看就是 4。这说明行列式值不依赖于“能不能除以主元”它天然就是一个整数再取模。于是我们需要换一套工具只允许用“第 3 条初等变换”和“第 1 条交换行”来完成消元因为这两条完全不需要求逆元。这就是所谓“不除法消元”的起点。2.3 辗转相除消元法的核心流程把问题局部化在消去第 i 列时我们面对的是第 i 行和第 j 行。当前这两行在第 i 列上的两个数分别是 x a[i][i] 和 y a[j][i]。我们的目标是在不除以任何数的前提下把 y 变成 0。这个过程本质上就是欧几里得算法。回想一下 gcd(x, y) 是怎么求的如果 x y交换然后用较大的减去若干倍较小的于是较大的变小。重复到最后一个变成 0另一个就是 gcd。矩阵里的“减法”对应的是第 3 条初等变换某行减另一行的倍数行列式不变“交换”对应第 1 条行列式变号。于是每一轮我们做三个动作计算 t x / y 的整数商对第 i 行执行 a[i] a[i] - t × a[j]交换第 i 行和第 j 行同时把答案 ans 取反。为什么交换后还要继续循环因为做完这轮后数对从 (x, y) 变成了 (y, x - t×y)而欧几里得算法正是通过不断翻转两个数来保证较大的数能迅速变小。重复直到 y 0此时第 j 行在第 i 列清零而第 i 行第 i 列的位置上留着一个非零的主元。整个过程里行列式变的只有符号。行倍加不改变行列式交换整行改变一次符号。我们把符号变化全部累计到 ans 变量里最后再把对角线乘积乘上去得到的就是原行列式乘上一个正负号也就是答案。我习惯用一句话记这个算法“用欧几里得的思路做高斯消元把除法全部换成加减和交换。”这也是 P7112 模板题的精髓。注意t 是整数除法的向下取整商。当 x y 时t 0此时“行倍加”这一步不会改任何数但“交换”照样会发生符号照样要翻。不少人第一次写就漏在这里结果样例都过不了。3. 完整 C 模板代码与逐行注释3.1 常规模板代码先把可直接提交的 C 模板贴出来我用的是 O(n^3 log mod) 的辗转相除写法n 最大 600 时在洛谷上可以稳定通过#include bits/stdc.h using namespace std; typedef long long ll; const int N 610; ll a[N][N]; int n; ll mod; ll det(int n, ll mod) { ll ans 1; for (int i 1; i n; i) { // 寻找当前列非零的行作为主元行 int pivot -1; for (int r i; r n; r) { if (a[r][i] % mod) { pivot r; break; } } if (pivot -1) return 0; // 整列全为 0行列式必为 0 // 若主元不在第 i 行交换过去同时符号取反 if (pivot ! i) { for (int j i; j n; j) swap(a[i][j], a[pivot][j]); ans (mod - ans) % mod; } // 用辗转相除消去第 i 行下面的所有 a[r][i] for (int r i 1; r n; r) { while (a[r][i] % mod) { ll t a[i][i] / a[r][i]; if (t) { for (int j i; j n; j) { a[i][j] (a[i][j] - t * a[r][j]) % mod; if (a[i][j] 0) a[i][j] mod; } } for (int j i; j n; j) swap(a[i][j], a[r][j]); ans (mod - ans) % mod; } } // 对角线主元乘进答案 ans ans * a[i][i] % mod; } return ans % mod; } int main() { scanf(%d%lld, n, mod); for (int i 1; i n; i) for (int j 1; j n; j) { scanf(%lld, a[i][j]); a[i][j] % mod; } printf(%lld\n, det(n, mod)); return 0; }3.2 逐段解释主元查找外层 for 循环每进入第 i 轮时前 i-1 列已经全部变成上三角结构接下来只需要处理第 i 到 n 列。先在 i 到 n 行里找第 i 列非零的行。如果找不到说明第 i 列全是 0行列式等于 0直接返回 0。这一步非常重要很多人忽略等到后面除零或者答案全员变 0 时才反应过来。主元交换把找到的非零行换到第 i 行。交换整行会让行列式变号所以 ans 要取反。注意我写的ans (mod - ans) % mod;是一种常见的取反写法等价于ans -ans;但能保证非负。如果你喜欢写ans -ans;也没问题最后结果再取一次模即可。辗转相除内层循环这是整个算法的核心。当第 r 行当前列还有值时反复执行用整数除法算出 t如果 t 非零做行倍减第 i 行减去 t 倍第 r 行无论如何都交换第 i 行和第 r 行ans 取反。为什么这里即使 t 0 也要交换因为 t 0 说明 a[i][i] 比 a[r][i] 小此时数对需要交换位置才能继续接下来的欧几里得过程。举个例子数对是 (2, 5)t 2/5 0不交换的话永远没法消交换后变成 (5, 2)t 5/2 2就能继续减了。这一步本质上是在模拟欧几里得的“大小翻转”。对角线乘入 ans第 i 轮结束后第 i 行在第 i 列的位置 a[i][i] 就是消元后的主元。把它乘进 ans。注意这里主元必须是非零的如果它在取模后变成 0说明这一列虽然处理完但最终主元与模数不互质导致行列式为 0。返回 0 也是对的。不过在上面的代码里如果 pivot 找到的非零值在后续消元中被模 p 变成了 0最后一步 ans ans * 0 0结果一样所以不用额外判断。4. 正确性验证与复杂度分析4.1 手工推演两个小规模案例光看算法容易觉得抽象我用手工推演的方式验证一遍。先看 n 2p 6矩阵[[2, 3], [4, 5]]行列式按定义算2×5 - 3×4 10 - 12 -2模 6 等于 4。下面用模板跑一遍。ans 1i 1找 pivot第 1 行第 1 列是 2非零直接用第 2 行消元a[2][1] 4非零进入 whilet 2 / 4 0行倍减不做交换第 1、2 行矩阵变成 [[4,5],[2,3]]ans 5再循环a[2][1] 2非零t 4 / 2 2第 1 行减 2 倍第 2 行[4 - 2×2, 5 - 2×3] [0, -1]取模 6 得 [0, 5]交换第 1、2 行矩阵变成 [[2,3],[0,5]]ans 1a[2][1] 0跳出 whileans × a[1][1] 2得 ans 2i 2没有需要消的ans × a[2][2] 5得 ans 10模 6 得 4。和手算一致。关键是符号整个过程经历两次交换负负得正所以行列式没有被翻转。再看一个模数为合数的 2×2 例子p 9矩阵 [[1,2],[3,4]]行列式 4 - 6 -2 ≡ 7。ans 1i 1a[1][1] 1对 r 2a[2][1] 3t 1 / 3 0交换行得到 [[3,4],[1,2]]ans 8a[2][1] 1t 3 / 1 3第 1 行减 3 倍第 2 行3 - 3×1 04 - 3×2 -2 ≡ 7于是第 1 行变成 [0,7]交换行矩阵变成 [[1,2],[0,7]]ans 1a[2][1] 0跳出ans × a[1][1] 1i 2ans × a[2][2] 7最终 ans 7。正确。这两个例子足够说明算法不仅能算符号处理也准确。4.2 复杂度与工程优化朴素版辗转相除消元的复杂度是 O(n^3 log mod)。原因是最外层 n 轮循环内层处理剩下 O(n) 行每次 while 迭代类似欧几里得取模最多 O(log mod) 次而行向量长度为 O(n)。n 600 时理论上限大约 2.16 亿次操作再乘以 log 因子听上去挺吓人。但实际运行远低于这个上限因为 while 循环的平均迭代次数远比 log mod 小而且一旦某行当前列变成 0 就立刻结束。洛谷上 n 600 的数据这个模板在 2 秒内跑完是没问题的。如果还是担心常数可以做几处优化行交换时从第 i 列开始不要从第 1 列开始。因为前 i-1 列已经是上三角结构交换它们只会浪费 O(n) 时间。少用%运算符。比如行倍减时先算出差值如果为负再补一次 mod而不是无脑(a mod) % mod。读入用scanf或快速读入不要用 cin 关同步。主元查找失败提前 return避免无意义计算。这些优化加起来实测能比朴素写法快 20% 到 40%。5. 实战踩坑WA、TLE、负数的三座大山5.1 最典型的五种错误我把这个模板反反复复交过很多次也在竞赛群里看过别人踩坑常见的错误基本集中在五个地方错误类型现象原因解决方案符号漏翻样例不过答案差正负号pivot 交换或辗转相除中漏写 ans 取反每次交换整行后立即ans -ans主元为 0 还在消运行时除零或答案全 0没有检查当前列是否全零先找 pivot找不到直接 return 0负取模输出负数行倍减后得到负数直接存下每次更新后调整到 [0, mod)使用逆元合数模数下答案错乱还在用质数模板的老思路改用欧几里得行变换不除主元复杂度爆炸TLEwhile 循环中每次都用% mod更新且从第 1 列开始优化循环边界并减少取模第五个问题尤其值得展开。很多人把内部循环写成for (int j 1; j n; j) a[i][j] (a[i][j] - t * a[r][j]) % mod;这从正确性上没错但多了很多无用操作而且当 t × a[r][j] 很大时C 取模要计算昂贵的除法。高常数叠加起来600 的数据很可能直接 TLE。正确的边界是从第 i 列开始因为前 i-1 列已经为 0做不做都不影响结果。5.2 对拍调试的傻瓜也能学会的方法模板代码写完以后我强烈建议写一个暴力求行列式的对拍程序。暴力版不需要快但一定要按定义写对ll detBrute(ll mat[][N], int n, ll mod) { if (n 1) return mat[1][1] % mod; ll res 0; for (int j 1; j n; j) { // 提取余子式 ll sub[N][N]; int pi 2, pj 1; for (int r 2; r n; r) { pj 1; for (int c 1; c n; c) { if (c j) continue; sub[pi-1][pj] mat[r][c]; pj; } pi; } ll term mat[1][j] * detBrute(sub, n - 1, mod) % mod; if (j % 2 1) res (res term) % mod; else res (res - term mod) % mod; } return res % mod; }这个暴力只适合 n ≤ 7 的小数据但配合随机矩阵生成器能在 10 分钟内帮你找出大部分符号和取模错误。对拍时随机生成 n 在 1 到 6 之间的矩阵模数随机选合数比如 4、6、8、9、12多跑几千组。一旦模板和暴力结果全对齐基本可以放心提交。6. 模板的复用矩阵树定理与生成树计数6.1 把行列式模板套到生成树计数P7112 这个模板最大的价值不是“会做这一道题”而是能直接迁移到矩阵树定理的场景。矩阵树定理说的是一个无向图的生成树数量等于它的拉普拉斯矩阵去掉任意一行一列后剩下矩阵的行列式。很多计数题最后都要求对某个模数取模而模数很可能不是质数。比如题目给你一张图让你输出生成树数量 mod 12345678这个数显然不是质数。如果你只会质数版本的高斯消元到了这一步就只能干瞪眼。复用方式很简单先读入图构造 n×n 的拉普拉斯矩阵删掉最后一行最后一列得到一个 (n-1)×(n-1) 的矩阵然后直接调用上面的 det 函数。// 假设已经建好 n 个点的拉普拉斯矩阵 lap[N][N] int m n - 1; ll ans det(m, mod); printf(%lld\n, ans);你甚至不需要修改 det 函数因为模板只依赖矩阵本身和模数不关心矩阵是怎么来的。这也是“模板题”名称的真正含义算法本身是独立的、可嵌入任何题目的。6.2 模板化的注意点复用时有几个细节要注意矩阵的编号需要手动调整。如果你读入的图是 1 到 n 编号构造拉普拉斯矩阵时可以直接把第 n 行第 n 列去掉剩下三维数组依然是 1 到 n-1。无向图每条边会让拉普拉斯矩阵的对角线 1且两个端点互相 -1。有重边要对应累加自环通常忽略具体看题目要求。如果图不连通行列式结果是 0det 函数里整列为 0 时直接 return 0正好命中。有向图生成树计数用的是有向拉普拉斯矩阵构造方式不同但底层行列式求值逻辑一模一样。这个模板还可以用在其他需要计算行列式的地方比如伴随矩阵、基尔霍夫矩阵、某些概率生成函数问题甚至线性代数里判断矩阵是否可逆。把这些应用场景串起来你就发现学一个模板题能辐射出至少五类进阶题。7. 写在模板后面的个人经验最后再分享一个自己踩出来的经验不要一上来就把这份模板背下来先亲手把 2×2 和 3×3 的手工推演各做一遍再抄代码最后再不看代码默写一遍。我当初就是先背代码结果换了一道模数特别小的题直接懵了回来重新理解了欧几里得消元的每一步才真正过关。模板只是结果真正值钱的是“为什么每次交换都要翻符号”和“为什么 t 会是整数除法”这两个问题的答案。把这些想透了P7112 这道题才算真正吃进肚子里以后遇到任何非质数模数的行列式你都能面不改色地掏出这套模板。