矩阵快速幂

一、暴力递推的死局：当 $n$ 暴涨到 $10^{18}$

递推是我们在入门算法时最熟悉的朋友。比如最经典的斐波那契数列：

让你求第 $n$ 项，闭着眼睛都能写出一个 for 循环：开两个变量滚动累加，一步一步往后推，$O(n)$ 轻松解决。

但如果出题人把数据范围改成了 $n \le 10^{18}$ 呢？

计算机每秒执行的基本运算大约是 $10^8$ 次。如果你老老实实循环 $10^{18}$ 次，程序需要跑几百年，评测机只会无情地甩给你一个超时（TLE）大礼包。

痛点根源：递推的本质是“单步更新”，每算一项就必须消耗一次循环。

破局思路：回想我们在初学整数求幂 $a^n$ 时，如果从 $1$ 连乘 $n$ 次也是 $O(n)$；但借助二进制拆分的快速幂，我们可以把 $a \times a$ 打包成 $a^2$，再打包成 $a^4, a^8\dots$，只用 $\log n$ 次乘法就搞定。

那能不能把递推中“求下一项”的操作，也打包成一次乘法？

只要能把状态转移写成乘法的形式，我们就能把快速幂套上去，将庞大的 $10^{18}$ 步转移压成几十轮矩阵运算！这个充当“打包工具箱”的数学武器，就是矩阵。

二、破局魔法：把“一次转移”打包成矩阵乘法

线性代数的理论很深，但作为信奥选手，我们先不要去背那些抽象的行列式与特征值。今天我们只从状态转移的物理视角来看矩阵。

1. 什么是状态行向量？

在信奥代码中，我们最习惯将当前的所有状态横排放在一个行向量里。

假设系统里有两个变量 $x$ 和 $y$，当前状态就是 $[x, y]$。

现在我们要执行一个操作：第一个量保持不变，第二个量变成两者的和。也就是：

我们把这个变化规则写成一个表格（矩阵），排成两行两列：

为什么是这个矩阵？我们来拆开看它的物理意义：

• 第一列负责拼出“新的第一个数”：旧的 $x \times 1 + y \times 0 = x$。
• 第二列负责拼出“新的第二个数”：旧的 $x \times 1 + y \times 1 = x + y$。

记忆法则（关键）： 矩阵中的每一行，代表旧状态里的某一个变量； 矩阵中的每一列，代表新状态里的某一个变量的制造配方！

我们把这个矩阵记为 $A$。再假设有另一个操作 $B$：“第一个量变成两者之和，第二个量保持不变”，它的矩阵就是：

如果从初始状态 $[2, 3]$ 出发，先做操作 $A$ 得到 $[2, 5]$，紧接着做操作 $B$ 得到 $[7, 5]$。这一连串过程可以写成：

2. 矩阵乘法的物理意义：中间人 $k$ 的搭桥

能不能把操作 $A$ 和操作 $B$ 提前合成一个整体操作 $C = AB$？

数学上的矩阵乘法定义为：若 $A$ 是 $r \times s$ 的矩阵，$B$ 是 $s \times t$ 的矩阵，则相乘得到的 $C = AB$ 是一个 $r \times t$ 的矩阵，其中第 $i$ 行第 $j$ 列的元素为：

看懂这个式子，就抓住了矩阵乘法的灵魂：

• $k$ 是中间人！从旧状态的第 $i$ 项，通过第一阶段转移到中间状态的第 $k$ 项，贡献是 $A_{ik}$；
• 接着从中转项 $k$ 通过第二阶段转移到最终的第 $j$ 项，贡献是 $B_{kj}$；
• 把所有可能的中间人 $k$ 串联起来相乘，再全部加在一起，就得到了从 $i$ 到 $j$ 的总贡献！

正因为必须有“中间人”来对接，所以前一个矩阵的列数必须等于后一个矩阵的行数，否则中转接口完全对不上。

我们亲手算一下刚才的 $A \times B$：

我们用合成后的矩阵来算一次：$[2, 3] \begin{pmatrix} 2 & 1 \\ 1 & 1 \end{pmatrix} = [2 \times 2 + 3 \times 1, 2 \times 1 + 3 \times 1] = [7, 5]$。结果完全一致！

如果倒过来算 $BA$，你会发现：

从 $[2, 3]$ 出发做 $BA$，结果是 $[5, 8]$。 致命警告：矩阵乘法不满足交换律（$AB \ne BA$）。先做谁后做谁，顺序绝对不能乱。

3. 为什么能用快速幂？结合律与单位矩阵

虽然矩阵乘法不能随意交换前后顺序，但它严格满足结合律：

只要顺序不变，先算哪一对完全自由！

如果我们要连续做 $n$ 次相同的状态转移 $T$，也就是求：

既然满足结合律，我们就可以把 $T \times T$ 合并为 $T^2$，把 $T^2 \times T^2$ 合并为 $T^4$……这就是经典的快速幂。

在整数快速幂里，乘法的初始基础是 $1$；而在矩阵世界里，充当这个“1”的是单位矩阵 $I$： 主对角线全为 $1$，其余所有位置全为 $0$。

单位矩阵的物理意义极其纯粹：什么都不做，让所有状态完美保持原样。任何矩阵乘上单位矩阵都等于自身（$TI = IT = T$），且规定 $T^0 = I$。

三、经典第一战：斐波那契数列的极速突围

1. 状态设计与系数推导

现在我们来攻克斐波那契数列：$F_{t+2} = F_{t+1} + F_t$。

想要算出下一项，光知道当前项 $F_t$ 是不够的，必须同时握有相邻的两项。因此，我们的灵魂状态向量必须打包两个变量：

我们的目标是找到一个 $2 \times 2$ 的转移矩阵 $T$，使得：

按照前面讲的“看列填系数”法则：

1. 新状态的第 1 格需要变成 $F_{t+1}$： 它等于 $0 \times F_t + 1 \times F_{t+1}$，所以矩阵第 1 列填 $\begin{pmatrix} 0 \\ 1 \end{pmatrix}$。
2. 新状态的第 2 格需要变成 $F_t + F_{t+1}$： 它等于 $1 \times F_t + 1 \times F_{t+1}$，所以矩阵第 2 列填 $\begin{pmatrix} 1 \\ 1 \end{pmatrix}$。

把两列拼起来，转移矩阵 $T$ 破茧而出：

我们手推演练一下从初态 $S_0 = [F_0, F_1] = [0, 1]$ 出发的变化过程：

• $[0, 1] T = [1, 1] = [F_1, F_2]$
• $[1, 1] T = [1, 2] = [F_2, F_3]$
• $[1, 2] T = [2, 3] = [F_3, F_4]$
• $[2, 3] T = [3, 5] = [F_4, F_5]$

转移一次下标推进一步，因此转移 $n$ 次后的状态就是：

2. 答案坐标怎么取？别靠猜！

设矩阵 $T^n$ 算出来的结果是：

那么最终的状态行向量为：

对比左右两边，立即得到：

• $F_n = c$，也就是 $T^n$ 中第 1 行第 0 列的元素（采用 0-based 数组下标即 R.a[1][0]）！
• $F_{n+1} = d$，也就是 R.a[1][1]。

考场上千万不要去死记最后答案是取 a[0][0] 还是 a[1][0]。花 5 秒钟写出初态行向量 $[0, 1]$ 与结果方阵的一行乘法，位置清晰明了，想错都难。

3. 核心模板一：斐波那契极速求解

输入输出协议与范围：

• 输入：一个非负整数 $n$（$0 \le n \le 10^{18}$）。
• 输出：$F_n \bmod 1000000007$ 的值。
• 样例测试：
  • 输入 10，输出 55
  • 输入 0，输出 0
• 复杂度：矩阵大小为常数 $2$，单次乘法运算仅 $2^3 = 8$ 次，时间复杂度严格为 $O(\log n)$，空间复杂度 $O(1)$。

#include<bits/stdc++.h>
#define int long long
using namespace std;
const int MOD=1000000007;

struct Mat{
	int a[2][2];
	Mat(){
		memset(a,0,sizeof(a));
	}
};

// 矩阵乘法 C = A * B
Mat mul(const Mat& A,const Mat& B){
	Mat C;
	for(int i=0;i<2;i++){
		for(int k=0;k<2;k++){
			if(!A.a[i][k]) continue; // 剪枝小优化
			for(int j=0;j<2;j++){
				C.a[i][j]=(C.a[i][j]+A.a[i][k]*B.a[k][j])%MOD;
			}
		}
	}
	return C;
}

// 矩阵快速幂 R = A^k
Mat power(Mat A,int k){
	Mat R;
	R.a[0][0]=R.a[1][1]=1; // 初始化为单位矩阵 I
	while(k){
		if(k&1) R=mul(R,A);
		A=mul(A,A);
		k>>=1;
	}
	return R;
}

void solve(){
	int n;
	if(!(cin>>n)) return;
	Mat T;
	// 构造转移矩阵 T
	T.a[0][0]=0; T.a[0][1]=1;
	T.a[1][0]=1; T.a[1][1]=1;
	
	Mat R=power(T,n);
	// 初始状态 [0, 1] 右乘 R，提取第 1 行第 0 列
	cout<<R.a[1][0]<<'\n';
}

signed main(){
	ios::sync_with_stdio(0),cin.tie(0);
	solve();
	return 0;
}

四、状态构造的降维拆解：变式递推如何搭积木？

学会了斐波那契，最怕换一道递推题就两眼抹黑。其实构造转移矩阵有一套通用的“搭积木流程”：

1. 确定状态维度：写出递推式，看要推出下一时刻，最少需要抓住历史上的哪几个量。
2. 写出行向量：把这些量排成一行 $[x_1, x_2, \dots]$。
3. 按列配平系数：针对每一个目标新变量，看它由旧变量的几倍相加而成，依次把系数填进对应的列。

1. 变式 1：带系数的二阶线性递推

递推式：$f_{t+2} = 2f_{t+1} + 3f_t$，给定初值 $f_0 = 1, f_1 = 2$。

我们依然维护当前两项 $[f_t, f_{t+1}]$，下一步目标是 $[f_{t+1}, 3f_t + 2f_{t+1}]$：

• 新第 1 格：$0 \times f_t + 1 \times f_{t+1}$ $\to$ 第 1 列填 $\begin{pmatrix} 0 \\ 1 \end{pmatrix}$
• 新第 2 格：$3 \times f_t + 2 \times f_{t+1}$ $\to$ 第 2 列填 $\begin{pmatrix} 3 \\ 2 \end{pmatrix}$

转移矩阵立得：

我们手算验算： 初态为 $[1, 2]$。 $[1, 2] \begin{pmatrix} 0 & 3 \\ 1 & 2 \end{pmatrix} = [2, 1 \times 3 + 2 \times 2] = [2, 7]$，对应 $f_1=2, f_2=7$。 $[2, 7] \begin{pmatrix} 0 & 3 \\ 1 & 2 \end{pmatrix} = [7, 2 \times 3 + 7 \times 2] = [7, 20]$，对应 $f_2=7, f_3=20$。 完全吻合！求第 $n$ 项时，最终答案为初态乘以 $T^n$：

此时 $f_n$ 就要取 $1 \times (T^n)_{0,0} + 2 \times (T^n)_{1,0}$，不能再只看单个格子了。

2. 变式 2：常数项骚扰 —— 引入“常数 1”扩维

递推式：$f_{t+1} = 2f_t + 3$。

多了一个带尾巴的常数 $+3$，怎么处理？ 扩维技巧：如果一个量在递推中不是变量，我们也可以把它当成一个“永远不变的特殊变量”塞进状态向量！ 我们构造状态向量为 $[f_t, 1]$：

• 下一步的新第一格是 $2f_t + 3 \times 1$；
• 下一步的新第二格常数必须维持为 $1$（即 $0 \times f_t + 1 \times 1$）。

矩阵立刻写出：

验算：$[f_t, 1] \begin{pmatrix} 2 & 0 \\ 3 & 1 \end{pmatrix} = [2f_t + 3, 1]$，完美契合！

3. 变式 3：多阶递推的“传送带”法则

递推式：$f_t = f_{t-1} + f_{t-2} + f_{t-3}$，初值为 $f_0=1, f_1=1, f_2=2$。

需要用到过去连续三项，状态向量就扩为三维 $[f_t, f_{t+1}, f_{t+2}]$。 下一步目标是 $[f_{t+1}, f_{t+2}, f_t + f_{t+1} + f_{t+2}]$：

• 前两格就像“传送带”一样往左平移一格：
  • 新第 1 格取旧第 2 格：列为 $\begin{pmatrix} 0 \\ 1 \\ 0 \end{pmatrix}$
  • 新第 2 格取旧第 3 格：列为 $\begin{pmatrix} 0 \\ 0 \\ 1 \end{pmatrix}$
• 新第 3 格把旧三项全部相加：列为 $\begin{pmatrix} 1 \\ 1 \\ 1 \end{pmatrix}$

组合起来就是优雅的伴随矩阵：

4. 递推前缀和：把“历史账本”塞进矩阵

很多时候，题目不仅要求我们算出第 $n$ 项，还会变本加厉地要求前 $n$ 项的和：$S_n = \sum_{i=0}^n F_i$。

如果单独算每一项再加起来，运算量又退化成了 $O(n)$。怎么办？ 扩维策略：把前缀和 $S_t$ 当成一个普通的状态变量，一起打包塞进矩阵里！

以斐波那契数列为例，我们需要维护的状态从两个扩展到了三个：$[F_t, F_{t+1}, S_t]$。 下一步的目标是推导出 $[F_{t+1}, F_{t+2}, S_{t+1}]$。 其中，新的前缀和 $S_{t+1} = S_t + F_{t+1}$。

按照“看列填系数”法则：

1. 新第 1 格 ($F_{t+1}$)：等于 $0 \times F_t + 1 \times F_{t+1} + 0 \times S_t$，所以第 1 列填 $\begin{pmatrix} 0 \\ 1 \\ 0 \end{pmatrix}$。
2. 新第 2 格 ($F_{t+2}$)：等于 $1 \times F_t + 1 \times F_{t+1} + 0 \times S_t$，所以第 2 列填 $\begin{pmatrix} 1 \\ 1 \\ 0 \end{pmatrix}$。
3. 新第 3 格 ($S_{t+1}$)：等于 $0 \times F_t + 1 \times F_{t+1} + 1 \times S_t$，所以第 3 列填 $\begin{pmatrix} 0 \\ 1 \\ 1 \end{pmatrix}$。

把这三列并排拼起来，得出的扩维转移矩阵为：

$T = \begin{pmatrix} 0 & 1 & 0 \\ 1 & 1 & 1 \\ 0 & 0 & 1 \end{pmatrix}$

我们从初态 $[F_0, F_1, S_0] = [0, 1, 0]$ 出发手算推演一步： $[0, 1, 0] \begin{pmatrix} 0 & 1 & 0 \\ 1 & 1 & 1 \\ 0 & 0 & 1 \end{pmatrix} = [1, 1, 1]$。 结果对应 $F_1=1, F_2=1, S_1=1$（$S_1 = F_0 + F_1 = 0 + 1 = 1$），完美契合！ 通过仅仅增加一个维度，我们不仅顺手求出了前缀和，还能继续用同一套快速幂代码在 $O(\log n)$ 内秒杀问题。

5. 周期性转移：将一组规律打包为“时间胶囊”

有时候，题目给的转移规律不是一成不变的，而是周期性变化的。 比如：第 1 秒按照矩阵 $A$ 转移，第 2 秒按 $B$ 转移，第 3 秒按 $C$ 转移，然后循环往复（$A, B, C, A, B, C \dots$）。 问经过 $k$ 秒后的状态。

这还能用快速幂吗？当然可以！虽然每个时刻矩阵不同，但矩阵乘法满足结合律，我们可以把一整个周期的转移合并打包成一个总矩阵：

$T_{period} = A \times B \times C$

行向量的物理直觉优势：这里必须强调一下我们坚持使用行向量的巨大优势！ 因为状态在左侧，操作在右侧，状态演化的书写顺序是 $[State] \times A \times B \times C$。从左往右读，时间流动的方向与阅读顺序完全一致！如果你用的是列向量，式子就得反过来写成 $C \times B \times A \times \begin{pmatrix} State \end{pmatrix}$，非常容易把顺序乘反导致全盘皆输。

有了周期总矩阵 $T_{period}$，对于 $k$ 秒的演化，我们只需要：

1. 算出完整的周期数 $p = k / 3$。
2. 对周期总矩阵做快速幂：$R = (T_{period})^p$。
3. 算出剩下的零头步数 $rem = k \bmod 3$。
4. 最终状态就是：$Initial\_State \times R \times (\text{如果 } rem \ge 1 \text{ 再乘 } A) \times (\text{如果 } rem \ge 2 \text{ 再乘 } B)$。

6. 边界警告：什么递推绝对不能直接打包？

这里的矩阵快速幂适合加速系数固定的线性递推，常数项可以像上一例一样，扩一维 $1$ 来吸收。

1. 非线性运算：如果递推式是 $f_t = f_{t-1} \times f_{t-2}$，变量之间发生了乘积，这无法用固定的系数做线性加权相加，矩阵直接失效（某些乘积递推可以通过取对数转化为加法，但那是另一回事）。
2. 动态变化的转移系数：如果递推系数每一步都随时间 $t$ 随机乱变，那每一步乘的都是完全不同的矩阵，无法使用同一底数的快速幂。

五、跨界融合：图论游走与邻接矩阵的连击魔法

现在，我们要把视线从数论递推跳跃到一个看起来截然不同的领域：图论。

场景：给出一张有向图，问：从起点 $s$ 出发，恰好走 $k$ 条边到达终点 $t$，一共有多少种不同的方案？ (注意：这里允许重复经过顶点和边，允许有重边和自环，这种在图上漫步的路径称为“游走”)

1. 物理意义的升华：中间点就是中转站

我们用邻接矩阵 $A$ 来存图：$A_{ij}$ 表示从节点 $i$ 到节点 $j$ 的有向边条数。

• 如果没有边，$A_{ij} = 0$；
• 如果有 1 条边，$A_{ij} = 1$；
• 如果有 2 条重边，$A_{ij} = 2$。

走 1 步的方案数，显然就是矩阵 $A$ 本身。

那走 2 步呢？ 从 $i$ 出发走 2 步到 $j$，必定要经过一个中间节点 $u$。

• 第一步从 $i$ 到 $u$ 有 $A_{iu}$ 条路；
• 第二步从 $u$ 到 $j$ 有 $A_{uj}$ 条路；
• 根据乘法原理，经过中转站 $u$ 的路径数是 $A_{iu} \times A_{uj}$。 枚举所有可能的中转节点 $u$，走 2 步的总方案数就是：

看看这个式子——这不就是矩阵乘法 $A \times A = A^2$ 的定义式吗！

在代数推导里抽象的“下标 $k$”，在图论中竟然瞬间具象化为了图上的中转节点 $u$！ 同理：

• 走 3 步就是 $A^3$；
• 恰好走 $k$ 步的方案数矩阵，就是邻接矩阵的 $k$ 次方 $A^k$！ 从起点 $s$ 到终点 $t$ 恰好走 $k$ 步的方案数，就是 $(A^k)_{st}$。

2. 三点小图的手推演练

给出一个具体的小图：

• 节点 1 到节点 2 有 2 条边
• 节点 2 到节点 1 有 1 条边
• 节点 2 到节点 3 有 1 条边
• 节点 3 到节点 3 有 1 条自环边

邻接矩阵 $A$ 写出来是：

我们来问：从 1 出发恰好走 4 步到 3，有几种走法？

我们先手算 $A^2$：

再算 $A^4 = A^2 \times A^2$。我们只关心第 1 行第 3 列： $(A^4)_{1,3} = \sum_u (A^2)_{1,u} (A^2)_{u,3} = 2 \times 2 + 0 \times 1 + 2 \times 1 = 4 + 0 + 2 = 6$ 种！

我们去图上肉眼验证这 6 种是哪来的：

• 路径类型 A：$1 \to 2 \to 3 \to 3 \to 3$。前两步到 3，后面绕自环 2 次。因为第 1 步有 2 条平行边可选，所以有 $2 \times 1 \times 1 \times 1 = 2$ 种。
• 路径类型 B：$1 \to 2 \to 1 \to 2 \to 3$。第 1 步有 2 条边可选，第 3 步也有 2 条边可选，共有 $2 \times 1 \times 2 \times 1 = 4$ 种。 合计整整 6 种！矩阵乘法神不知鬼不觉地帮我们枚举了所有的重边与自环组合。

3. 走 0 步的深刻内涵：为什么对角线必须是 1？

当 $k = 0$ 时，从起点出发走 0 步：

• 如果起点就是终点（$s = t$），“原地发呆”本身就是一种合法的方案，答案是 $1$；
• 如果起点不同于终点（$s \ne t$），走 0 步绝不可能到达其他点，方案是 $0$。

这恰恰与单位矩阵 $I$ 的定义（主对角线为 1，其余为 0）在物理层面上完美统一！

4. 核心模板二：恰好走 $k$ 步的路径方案数

输入输出协议与范围：

• 第一行输入 5 个整数：$n, m, k, s, t$
  • 点数 $1 \le n \le 50$
  • 边数 $0 \le m \le 10^4$
  • 步数 $0 \le k \le 10^{18}$
  • 起点终点 $1 \le s, t \le n$
• 接下来 $m$ 行，每行两个整数 $u, v$，表示一条从 $u$ 到 $v$ 的有向边（可能存在自环与重边）。
• 输出一个整数：从 $s$ 到 $t$ 恰好走 $k$ 步的方案数对 $1000000007$ 取模的结果。
• 复杂度：点数 $n \le 50$，矩阵快速幂耗时 $O(n^3 \log k)$，建图耗时 $O(m)$，空间复杂度 $O(n^2)$。别只看指数小，矩阵维度也要一起估算。

样例输入：

3 5 4 1 3
1 2
1 2
2 1
2 3
3 3

样例输出：

6

#include<bits/stdc++.h>
#define int long long
using namespace std;
const int MOD=1000000007;
const int N=55;

int n,m,k,s,t;

struct Mat{
	int a[N][N];
	Mat(){
		memset(a,0,sizeof(a));
	}
};

// 矩阵乘法，只循环到有效点数 n
Mat mul(const Mat& A,const Mat& B){
	Mat C;
	for(int i=1;i<=n;i++){
		for(int k_mid=1;k_mid<=n;k_mid++){
			if(!A.a[i][k_mid]) continue; // 核心剪枝：无前驱连边直接跳过
			for(int j=1;j<=n;j++){
				C.a[i][j]=(C.a[i][j]+A.a[i][k_mid]*B.a[k_mid][j])%MOD;
			}
		}
	}
	return C;
}

// 矩阵快速幂
Mat power(Mat A,int p){
	Mat R;
	for(int i=1;i<=n;i++) R.a[i][i]=1; // 初始化为 n 阶单位矩阵
	while(p){
		if(p&1) R=mul(R,A);
		A=mul(A,A);
		p>>=1;
	}
	return R;
}

void solve(){
	if(!(cin>>n>>m>>k>>s>>t)) return;
	Mat A;
	for(int i=1;i<=m;i++){
		int u,v;
		cin>>u>>v;
		A.a[u][v]=(A.a[u][v]+1)%MOD; // 支持重边累加
	}
	Mat R=power(A,k);
	cout<<R.a[s][t]<<'\n'; // 输出从 s 到 t 恰好走 k 步的方案数
}

signed main(){
	ios::sync_with_stdio(0),cin.tie(0);
	solve();
	return 0;
}

六、避坑指南：取模安全与防溢出铁律

1. 为什么 long long 不会在累加时爆掉？

在上面的乘法模板中，内层循环的核心语句是：

C.a[i][j] = (C.a[i][j] + A.a[i][k] * B.a[k][j]) % MOD;

很多同学会担心：乘积会很大，再加起来会不会爆 long long？ 我们来看极值估算：

• 矩阵内每一个元素经过取模后，最大是 $MOD - 1 \approx 10^9$；
• 两个数相乘，最大约为 $(10^9)^2 = 10^{18}$；
• 加上已经取模过的 C.a[i][j]，总和最大约为 $10^{18} + 10^9$；
• 而有符号 64 位整型 long long 的上限是 $2^{63}-1 \approx 9.22 \times 10^{18}$！

因为我们每加一次就取一次模，运算过程牢牢被压在 $10^{18}$ 左右，绝对安全。 但是！如果你贪图常数小，把取模移到最外层（把内层所有 $n$ 项全部乘完加在一起再取模）： 当 $n = 50$ 时，$50 \times 10^{18} = 5 \times 10^{19} > 9.22 \times 10^{18}$，会瞬间爆 long long 导致溢出变成负数！考场上写这种“危险优化”极易爆零。

2. 模数的本质：矩阵快速幂需要质数吗？

很多同学做多了数论题，一看到取模就想费马小定理，以为 $MOD$ 必须是质数。 大错特错！ 矩阵快速幂在整个运行过程中，只涉及乘法和加法，根本没有除法（逆元）运算。根据同余理论，加法和乘法对任意正整数模数都完美保持同余性质。哪怕模数是一个合数（比如 $10^9$ 或者 $65536$），算法依然百分之百正确。

七、图论进阶：Min-Plus 乘法与恰好 K 条边最短路（选学）

(先修要求：掌握《最短路与状态建图》中的 Floyd 全源最短路思想)

前面我们利用矩阵乘法 $C_{ij} = \sum_k A_{ik} B_{kj}$ 求出了恰好走 $K$ 步的“路径方案数”。 如果出题人换个问法：从起点 $s$ 走恰好 $K$ 条边到终点 $t$，最短距离是多少？

我们重新审视传统乘法里的操作物理意义：

• $A_{ik} \times B_{kj}$：代表两段路径组合方案数的累乘；
• $\sum_k$：代表把经过所有可能中转站 $k$ 的方案数累加。

在最短路问题中，组合两段路径的逻辑变了：

• 经过中转站 $k$ 时，两段路径的长度是相加（$A_{ik} + B_{kj}$）；
• 在所有可能的中转站 $k$ 中，我们要挑一个总长度最小的（$\min_k$）。

我们大胆地把矩阵乘法的加法与乘法内核替换掉：

这就是鼎鼎大名的 Min-Plus（最小-加）矩阵乘法！ 最神奇的是，因为 $\min(a+b, a+c) = a + \min(b, c)$（数学上的分配律依然成立），所以这种魔改后的乘法完美保留了结合律！既然满足结合律，快速幂就能毫无保留地套上去，只需 $O(\log K)$ 次矩阵乘法，就能算出走 $K$ 条边的最短路；用朴素矩阵乘法时，总时间是 $O(n^3\log K)$。

1. 核心实现要点与避坑

使用 Min-Plus 乘法时，有几个致命的细节必须改掉：

1. 初始化：普通乘法找方案数，没边的格子是 0；Min-Plus 找最短路，没边的格子必须初始化为正无穷 INF。
2. 单位矩阵：普通乘法的单位矩阵是对角线为 1，其余为 0；Min-Plus 乘法中，单位矩阵代表“走 0 步”。原地不动距离是 0，而0步飞到别的点是不可达 INF。所以它的单位矩阵是对角线全为 0，其余全为 INF！
3. 防溢出：既然有加法，如果用系统的 INT_MAX，两个 INF 相加极易爆成负数。必须定义一个安全的无穷大，如 1e18。

2. 完整程序：恰好 K 条边最短路

自拟实战演练：

• 输入：第一行 $n, m, k, s, t$（$1\le n\le 50, 0\le m\le 1000, 0\le k\le 10^9, 1\le s,t\le n$），接下来 $m$ 行每行输入 $u, v, w$，表示 $u \to v$ 权重为 $w$ 的有向边（支持重边）。约定 $|w|\le10^6$，因此目标游走长度的绝对值至多 $10^{15}$，小于模板的 INF=1e18。
• 输出：恰好走 $k$ 条边从 $s$ 到 $t$ 的最短路。如果不可达输出 Impossible。
• 样例输入：
  text复制代码

  3 3 2 1 3
  1 2 5
  2 3 4
  1 3 1
  

• 样例输出：9（解释：1->2->3 长度 9 恰好是两条边；1->3 长度 1 虽然更短，但只有一条边，不符合“恰好两步”的要求）。

#include<bits/stdc++.h>
#define int long long
using namespace std;
const int INF = 1e18; // 安全的无穷大，防止 INF + INF 溢出
const int N = 55;

int n, m, k, s, t;

struct Mat {
	int a[N][N];
	Mat() {
		// Min-Plus 初始化必须全是 INF
		for(int i=0; i<N; i++) {
			for(int j=0; j<N; j++) {
				a[i][j] = INF;
			}
		}
	}
};

// 魔改版：Min-Plus 矩阵乘法
Mat mul(const Mat& A, const Mat& B) {
	Mat C;
	for(int i=1; i<=n; i++) {
		for(int k_mid=1; k_mid<=n; k_mid++) {
			if(A.a[i][k_mid] == INF) continue; // 无路可走，直接跳过
			for(int j=1; j<=n; j++) {
				if(B.a[k_mid][j] == INF) continue;
				C.a[i][j] = min(C.a[i][j], A.a[i][k_mid] + B.a[k_mid][j]);
			}
		}
	}
	return C;
}

// Min-Plus 矩阵快速幂
Mat power(Mat A, int p) {
	Mat R;
	// Min-Plus 的单位矩阵：对角线为 0，其余为 INF
	for(int i=1; i<=n; i++) R.a[i][i] = 0; 
	
	while(p) {
		if(p & 1) R = mul(R, A);
		A = mul(A, A);
		p >>= 1;
	}
	return R;
}

void solve() {
	if(!(cin >> n >> m >> k >> s >> t)) return;
	Mat A;
	for(int i=1; i<=m; i++) {
		int u, v, w;
		cin >> u >> v >> w;
		A.a[u][v] = min(A.a[u][v], w); // 图上可能有重边，贪心保留最短的那条
	}
	
	Mat R = power(A, k);
	if(R.a[s][t] >= INF) cout << "Impossible\n";
	else cout << R.a[s][t] << '\n';
}

signed main() {
	ios::sync_with_stdio(0), cin.tie(0);
	solve();
	return 0;
}

八、渐进式实战练习题单

不要把矩阵快速幂当成冰冷的背诵模板，请通过以下阶梯练习建立肌肉记忆：

1. 阶段一：手动构造型练习

1. 变式二阶线性递推实测
   • 题目要求：输入一个整数 $n$（$0 \le n \le 10^{18}$），输出递推式 $f_0=1, f_1=2, f_{t+2}=2f_{t+1}+3f_t$ 的 $f_n \bmod 1000000007$。
   • 验证指引：输入 3 输出 20；输入 0 输出 1。
   • 突破点：修改第一份程序的转移矩阵系数，并注意当提取最终答案时，必须是初态向量 $[1, 2]$ 与幂矩阵的线性组合，不再是只取单个格子。

2. 阶段二：边界特判与图论机制

2. 零步与不可达特判测试

   • 题目要求：使用第二份图论程序，输入 2 0 0 1 1（没有边，步数为 0，起点终点均为 1），应当输出 1；将最后终点改成 2（输入 2 0 0 1 2），应当输出 0。
   • 突破点：如果你的代码两个都输出了 0，说明单位矩阵没初始化好；如果两个都输出了 1，说明忘记判断起点与终点是否重合。

3. 重边累乘爆炸测试

   • 题目要求：输入第一行 1 2 3 1 1，紧接着两行边均为 1 1（1 个点，2 条自环重边，走 3 步，从 1 到 1）。
   • 验证指引：输出应为 8（因为每一步都有 2 条边可选，走 3 步就是 $2^3 = 8$）。检验邻接矩阵是否正确累加了重边。

3. 阶段三：对拍与终极检验

4. 编写小范围暴力验证器
   • 训练要求：对于图论游走，当 $k \le 100$ 时，可以用最纯粹的动态规划做验证： 设 dp[step][u] 表示走 step 步到达节点 u 的方案数，每一层枚举每条边转移到 step + 1。
   • 突破点：拿随机小图将这个小暴力程序与你的矩阵快速幂程序进行对拍，重点检测带有重边、自环、步数为 0 以及死胡同不可达等极限情况。亲眼见证 $O(k)$ 的慢速递推与 $O(\log k)$ 的矩阵暴击给出完全一致的答案，你对“状态打包”的理解才算真正大功告成！