矩阵快速幂模板详解：从快速幂原理到代码实现

1. 题目理解与核心原理

1.1 这道模板题到底在考什么

P3390 这道题，几乎每个刷信奥题的人都绕不过去。题目名称里写着“模板”两个字，说明它就是拿来当标准件用的，理解了这道题，后续很多算法题都能直接套同一套思路。题目要求很直接：给定一个 n 阶方阵和一个正整数 k，求这个矩阵的 k 次幂，并对 1e9+7 取模。

初次接触的人可能会想：矩阵乘法不就是三层循环吗？直接乘 k 次不就行了？乍一看确实可以，但仔细算一下复杂度就傻眼了。单次矩阵乘法是 O(n^3)，乘 k 次就是 O(n^3 * k)。如果 n 取 100，k 取 1e9，这个计算量大概是 1e15 级别的操作，哪怕机器再快也得跑到地老天荒。所以这道题考察的核心，就是怎么把“连乘 k 次”变成“只做大约 log k 次矩阵乘法”，而这就是快速幂的思想。

把数字换成矩阵，把数字乘法换成矩阵乘法，快速幂的骨架完全没有变化。这也就是它被称为“模板”的原因：同一个二进制拆分的套路，既可以拿来算 2 的 100 次方，也可以拿来算某个递推矩阵的 10 亿次方。理解这一层，比单纯背下代码重要得多。

1.2 从普通快速幂到矩阵快速幂，关键一步是什么

我们先回忆一下普通快速幂。要计算 a^b，不是把 a 连乘 b 次，而是把 b 拆成二进制。比如 b = 13，二进制是 1101，对应 13 = 8 + 4 + 1，所以 a^13 = a^8 * a^4 * a^1。我们只需要从 a^1 开始，每次都把当前数平方，得到 a^2、a^4、a^8，根据二进制位决定要不要乘进答案。这样乘法次数从 b 次降到了约 2 * log2(b) 次。

矩阵快速幂完全同理。唯一的差别是：底数从数字 a 变成了矩阵 A，乘法从数字乘法变成了矩阵乘法，初始答案从 1 变成了单位矩阵 E。单位矩阵就是对角线全 1、其余全 0 的矩阵，任何矩阵乘上单位矩阵都等于它自己，相当于数字世界里的“1”。这一点必须想明白，因为很多人第一次写的时候习惯把结果矩阵初始化为全 0，最后发现算出来全是 0，这就是把数字快速幂里 res = 1 的初始化逻辑带偏了。

另一个容易忽略的点是矩阵乘法不满足交换律。A * B 和 B * A 通常不一样，所以代码里要特别注意累乘的顺序。不过在本题中，因为底数始终是同一个矩阵 A，而 A 的不同幂之间是可以交换的（A^p 乘 A^q 等于 A^q 乘 A^p），所以顺序不会影响最终结果。但如果以后把模板改造成带初始向量的场景，顺序就必须严格对齐转移矩阵和向量的左右位置。

1.3 复杂度与数据范围的预判

n 和 k 的范围决定了很多实现上的细节。n 最大可以到 100，k 最大可以到 1e9 这个量级。矩阵乘法单次是 O(n^3)，也就是最多 100 的三次方，100 万次乘加运算，而快速幂最多执行约 60 次平方和约 60 次条件乘法，整体下来大概是 1 亿级别的基础运算，C++ 在 1 秒内完全能跑完。

这给我们的启示是：这道题不追求极端复杂的常数优化，但要求每一步都不能做无谓的浪费。比如如果没写快速幂，复杂度崩盘；如果乘法函数里频繁动态分配内存，常数也会被拖慢。另外，乘法过程中涉及两个 1e9 级别的数相乘，乘积会达到 1e18，超过 int 的范围，所以矩阵元素必须用 long long 存。这个细节几乎每个初学者都会踩。

2. 核心实现与边界处理

2.1 矩阵存储方式：固定数组是信奥场的稳妥选择

实现矩阵快速幂，第一个选择就是用什么结构存矩阵。常见方案有二维数组、vector 套 vector、以及一维数组模拟二维。我在实际写题时反复对比过，信奥场景下最稳妥的就是固定大小的二维数组，比如long long mat[105][105]。

原因有两个。第一，编译期大小固定，内存分配在栈上或静态区，没有运行期动态分配的损耗和风险；第二，vector 套 vector 虽然写起来灵活，但每次乘法都要传递引用、访问时多次间接寻址，常数更大。n 最大只有 100，固定开 105 的数组绰绰有余，内存也就 100 * 100 * 8 字节，约 80KB，完全不用担心。

用固定数组需要注意函数传参的写法。二维数组传进函数时会退化成指向数组的指针，所以形参要写成ll mat[][105]或者ll (*mat)[105]，不能漏掉第二维的大小。我见过有人写成ll mat[][]，编译都过不了。如果觉得这种写法难看，也可以用一维数组加偏移量手动索引，但可读性会下降。对模板题来说，清晰比炫技更重要。

2.2 乘法函数里的一个隐蔽大坑：源与目标重叠

矩阵乘法本身不难，就是三层循环：C[i][j] = sum(A[i][t] * B[t][j])。但写成函数时，如果直接往 C 数组里累加，而 C 和 A 或 B 恰好是同一个数组，就会出问题。比如在快速幂里，经常会写mul(A, A, A)，意思是将 A 平方，如果乘法函数边读边写同一个 A，后面计算用的 A[i][t] 可能已经被覆盖了，结果全部错乱。

我常用的解决办法是在乘法函数内部定义一个临时数组 tmp，把计算结果先全部写进 tmp，循环结束后再整体复制回目标数组 C。这样无论调用方传进来的 A、B、C 是否指向同一个数组，都不会出现覆盖问题。代码看起来多了一个复制环节，但 n 只有 100，复制 1 万个 long long 的开销微乎其微，换来的是逻辑上的绝对安全。

这里还有一个经验：矩阵乘法三层循环的遍历顺序，理论上对缓存命中率有影响。i-j-t 的顺序因为 C[i][j] 是连续累加目标，A[i][t] 按行读，B[t][j] 按列读，通常表现不差。更激进的优化可以改成 i-t-j 顺序，让 B[t][j] 也按行读，但本题 n 很小，差别并不明显，先把正确性保住再谈优化。

2.3 快速幂主循环与单位矩阵初始化

快速幂的骨架我写在主函数里，逻辑和数字快速幂一模一样：

while (exp > 0) { if (exp & 1) matMul(res, res, A); matMul(A, A, A); exp >>= 1; }

exp & 1判断当前二进制最低位是否为 1，是就把当前底数矩阵乘进结果；然后不管最低位是几，都让底数自乘平方；最后把指数右移一位。举例来说，k = 5 的二进制是 101，第一次循环最低位是 1，把 A 乘进 res，然后 A 平方变成 A^2；第二次循环最低位是 0，不乘，A 平方变成 A^4；第三次循环最低位是 1，把 A^4 乘进 res，最终 res = A * A^4 = A^5。

单位矩阵因为对角线是 1，其余是 0，所以初始化可以直接写：

memset(res, 0, sizeof(res)); for (int i = 0; i < n; i++) res[i][i] = 1;

这里强调一下，如果 k 等于 0，while 循环一次都不执行，res 保持单位矩阵，直接输出就是正确答案，没有任何额外判断。很多模板题的数据会带 k = 0 的边界，这个初始化习惯能帮你自动躲过去。

2.4 取模运算和 long long 的细节

10^9+7 这个模数大家都很熟悉，但真正写起来有几个细节要留意。矩阵元素读进来时，虽然题目说元素本身可能在模数范围内，但为了绝对安全，读入后立刻取一次模。乘法累加时，tmp[i][j] = (tmp[i][j] + A[i][t] * B[t][j]) % MOD，这个写法能保证每一步都不会溢出 long long。

为什么中间结果不会溢出？因为 A[i][t] 和 B[t][j] 都小于 1e9+7，乘积小于 1e18，而 long long 能表示的范围大约是 9.22e18，所以单个乘积没问题。但注意，如果 tmp[i][j] 没有及时取模，累加 100 次后可能达到 1e20，那就炸了。所以必须每加一次就取模，或者累计到一定次数再取模。模板题直接每步取模最省心。

输出时也不要忘记每个元素都已经是对 MOD 取模的结果。真正比赛时经常有人因为读入时没取模、中间乘法泄漏溢出、输出前忘取模，三个小问题叠加，导致答案对不上。取模这件事宁可多写，不要少写。

3. 完整模板与优化实战

3.1 一个可以直接提交的完整模板

我把整套代码整理成一个标准模板，注释写清楚关键位置，方便直接抄去修改：

#include <bits/stdc++.h> using namespace std; typedef long long ll; const ll MOD = 1000000007LL; const int MAXN = 105; int n; ll k; ll base[MAXN][MAXN]; ll res[MAXN][MAXN]; void matMul(ll C[MAXN][MAXN], ll A[MAXN][MAXN], ll B[MAXN][MAXN]) { ll tmp[MAXN][MAXN] = {0}; for (int i = 0; i < n; i++) { for (int j = 0; j < n; j++) { for (int t = 0; t < n; t++) { tmp[i][j] = (tmp[i][j] + A[i][t] * B[t][j]) % MOD; } } } for (int i = 0; i < n; i++) { for (int j = 0; j < n; j++) { C[i][j] = tmp[i][j]; } } } void fastPow(ll A[MAXN][MAXN], ll exp) { memset(res, 0, sizeof(res)); for (int i = 0; i < n; i++) { res[i][i] = 1; } while (exp > 0) { if (exp & 1) { matMul(res, res, A); } matMul(A, A, A); exp >>= 1; } } int main() { scanf("%d%lld", &n, &k); for (int i = 0; i < n; i++) { for (int j = 0; j < n; j++) { scanf("%lld", &base[i][j]); base[i][j] %= MOD; } } fastPow(base, k); for (int i = 0; i < n; i++) { for (int j = 0; j < n; j++) { printf("%lld%c", res[i][j], j == n - 1 ? '\n' : ' '); } } return 0; }

这里有几个值得说明的设计。全局数组 base 和 res 直接定义在 main 外面，避免栈空间被大数组吃紧，也方便函数直接使用 n 变量，不用每次传参都带尺寸。matMul 内部用局部 tmp 数组，虽然占栈空间约 80KB，但函数返回就释放，没问题。

3.2 常数优化与代码可读性的平衡

不少选手喜欢在这里做优化，比如把% MOD改成if (tmp >= MOD) tmp -= MOD，因为模运算在 CPU 上比加减法慢不少。官方题解也常有人讨论这个问题。但我的观点是：对这道模板题，n 只有 100，每步取模的写法完全够过题，没必要为了微弱的常数提升牺牲正确性。真正的性能瓶颈如果存在，也多半在别处。

如果确实想优化，可以这样改：内层循环累加时先不取模，等内层 t 循环结束后统一取模。因为每个 tmp[i][j] 最多累加 n 次，每次乘积小于 1e18，乘 n 次也才 1e20，还是超了，所以这种优化只适用于乘积更小的场景。更稳妥的折衷方案是每累加 8 次取一次模，或者用 unsigned long long 的天然溢出配合取模，但这些都属于经验丰富之后再去玩的花活，新手阶段先把标准写法吃透。

可读性方面，我建议把 matMul 的形参写成 A、B、C，并在注释里标明“C = A * B”，这样以后改造成别的矩阵题时，一眼就能看出调用关系。我在实际刷题中吃过亏：把一个函数的参数顺序写成 A、C、B，结果调用的地方全乱套，排查了半天才发现是参数名误导了自己。

3.3 实测中的数据表现与性能体会

我在本机用随机生成的 100 阶矩阵测试过这个模板，k 取到接近 1e9 时，整个程序跑完大约在 0.2 秒到 0.4 秒之间，远小于一般题目的时间限制。这说明模板的常数很安全，不需要再做额外的循环展开或 SIMD 优化。

不过有另一个细节影响实际耗时：矩阵乘法的调用次数。快速幂中，除了平方运算，条件成立时还会多做一次 res 与 A 的乘法。如果二进制里 1 的位数很多，额外乘法次数就多。比如 k 接近 2^30-1 时，二进制里 1 的个数约为 30 个，总乘法次数约 60 次；如果 k 只有 2^30，二进制只有最高位一个 1，额外乘法只有一次。这个差异在 n 很小时无所谓，但 n 达到 500 以上的大矩阵题时就能感受到了。

4. 常见错误与调试实录

4.1 高频错误速查表

把容易踩的坑列成一个表，方便以后写题时快速对照：

这张表里的前三条我基本每次带新人都会讲一遍，因为它们覆盖了 80% 的初学者报错场景。尤其是单位矩阵初始化，很多人会拿数字快速幂的模板去套，数字里初始化 1，矩阵里应该初始化单位矩阵，结果不假思索地全填 0，最后输出全 0 才恍然大悟。

4.2 用手算小矩阵验证思路

调试这类算法题，最有效的办法不是盯着屏幕读代码，而是构造一个 2 阶小矩阵，用笔手算出幂次结果，再和程序输出对比。比如取矩阵 [[1, 2], [3, 4]]，先求它的平方，手算得到 [[7, 10], [15, 22]]，拿程序跑一遍看是否一致。如果这一步是对的，说明矩阵乘法函数本身没问题。

接着可以验证快速幂的逻辑。取 k = 3，手算 A^3 = A^2 * A，得到结果和程序对比。再取 k = 4，验证平方再平方的路径。这样从小规模开始逐层验证，比直接跑大数据猜测错误位置快得多。我通常还会额外测一个 k = 0，确认输出单位矩阵；再测一个 n = 1 的退化情形，也就是只有一个元素的矩阵，此时矩阵乘法实际上退化成数字乘法，快速幂应该直接输出底数的 k 次方取模。

4.3 对拍测试与随机数据

当程序能通过手算样例后，再进一步就是对拍。对拍的意思是写一个朴素的暴力程序，直接连乘 k 次得到结果，然后生成多组小范围的随机数据，同时跑暴力版本和快速幂版本，对比输出是否完全一致。小范围的 k 可以取 0 到 20，n 可以取 1 到 5，这样暴力乘法的次数很少，可以快速穷举很多组数据。

对拍脚本的思路很简单：生成随机 n、k 和矩阵元素，分别喂给两个程序，用 diff 比较输出文件。如果随机测了几百组都一致，基本可以确认算法核心没有逻辑错误。再回到题目的大数据范围测性能和最终答案。这一步看着麻烦，但每次都能在正式提交前抓到意想不到的问题，尤其是那种只在特定 k 值下才触发的边界 bug。

5. 扩展应用与刷题实战建议

5.1 从模板题到斐波那契与线性递推

矩阵快速幂最大的价值体现在递推加速。拿斐波那契数列来说，普通地推下去是 O(n)，但把它写成矩阵形式：

F(n+1) = F(n) + F(n-1)

可以用一个 2 阶转移矩阵表示从 [F(n), F(n-1)] 到 [F(n+1), F(n)] 的变化。那么求第 n 项，本质上就是求这个转移矩阵的 n 次幂。这样一来，O(n) 的递推就变成了 O(log n) 的矩阵快速幂。这个应用几乎是所有进阶题的基础，很多关于线性递推、DP 状态转移的题目，最后都能落到矩阵乘法上。

我刷过的一道模拟题是这样的：给定一个数列，每一项由前面三项线性组合得到，求第 k 项。直接递推会超时，但构造一个 3 阶转移矩阵，然后跑快速幂，轻松解决。做这类题有个通用套路：先把状态向量定下来，再根据递推关系写出转移矩阵，最后套模板求幂，乘上初始向量。

5.2 图论中的路径计数问题

另一个经典扩展是邻接矩阵的幂。一张图的邻接矩阵 A 中，A[i][j] 表示节点 i 到节点 j 是否有直接边。那么 A^k 中的 (i, j) 元素，恰好表示从 i 到 j 长度为 k 的路径数量加权的某种含义。如果边有边权，转移矩阵还可以换成加权矩阵。这背后是矩阵乘法的本质：C[i][j] = sum(A[i][t] * B[t][j])，枚举了中间节点 t，把左右两段路径拼接起来。

这类题目通常会把路径计数和取模结合，正好就是 P3390 模板的变体。我的建议是把模板题的代码封装成可复用的函数，比如matMul和fastPow保持通用，以后遇到路径计数、DP 加速、递推加速，直接调整矩阵的构造方式，核心函数一行都不用改。

5.3 刷题顺序与能力进阶路线

如果刚开始接触矩阵快速幂，我建议按这个顺序练手：先确保能不看参考写出一份完整模板，再拿斐波那契类的题目验证理解，接着尝试用矩阵优化一个已知的 DP，最后回到图论路径计数题。每一步都要求自己手推转移矩阵，而不是直接抄别人代码。

等你对这个套路足够熟练，可以进一步接触更高阶的内容：比如矩阵加速带有限制条件的递推、多组转移矩阵的乘积、以及结合线段树支持修改的进阶题。但无论题目多复杂，底层都是同一个快速幂骨架。所以我一直觉得，P3390 这道模板题值得反复重写几次，直到闭着眼睛都能把单位矩阵初始化、二进制拆分循环、tmp 缓冲这三件事写对为止。

个人经验与最后一点建议

矩阵快速幂是我刷题路上第一个真正感受到“算法改变计算量级”的模板。以前觉得算 10 亿次幂是天方夜谭，后来发现只需要几十次乘法。这个认知转变让我对后续很多优化算法都建立了信心。如果你也正在做刷题打卡，遇到这道模板题时别急着背代码，先自己推一遍二进制拆分的原理，手算一个 2 阶矩阵验证每一步，再动手实现。踩过的坑越多，印象越深，以后遇到递推加速、路径计数等进阶题时，你会感谢现在认真扎实的自己。