矩阵快速幂加速
递推→矩阵幂·O(k³log n)
本课摘要
矩阵快速幂加速课程回答“线性递推怎样写成矩阵并用快速幂加速”。内容以递推→矩阵幂·O(k³log n)为主线，配合逐步推导、可编辑演示、例题与练习形成可复查的学习闭环。
- 判断矩阵快速幂加速的适用条件与状态边界
- 围绕“递推→矩阵幂·O(k³log n)”推导转移与计算顺序
- 用演示、复杂度分析和配套题目校验实现
本课目录 · 0 节
正在整理目录…
当递推项数大到没法一格一格填
回到最朴素的斐波那契：，从 起逐项往上加。 只要 不大，一个 的循环就够——这正是 B 部分计数 DP 里数楼梯那一套。 可一旦题目把 抬到 （例题 P1962 就是），逐项递推要走近 九百亿亿 步，任何机器都算不完。
瓶颈很清楚：递推一步只跨一项，代价被死死锁在 。要提速，就得想办法一次跨过很多项。 突破口在于——斐波那契这类递推是线性的（新项是旧几项的线性组合，没有平方、没有取最值）。线性变换恰好可以写成矩阵乘法，而「重复施加同一个线性变换 次」就是求矩阵的 次幂—— 幂运算有快速幂（二进制倍增），能把 次压成 次。这一节就把这条「递推 → 矩阵 → 快速幂」的加速链讲透。
把递推写成矩阵乘法
关键一步：把「当前需要记住的几项」打包成一个状态向量。斐波那契的新项只用到前两项，于是取行向量 作状态。 我们要找一个固定的矩阵 ，让它右乘这个向量后，正好整体前移一步，得到 ：
怎么定出 的每个数？逐个输出分量看它由旧分量怎样线性组合：新向量第一个分量要等于 ， 即旧的两个分量各取一份相加 → 对应 的第一列是 ；新向量第二个分量要等于 ， 即只取旧的第一个分量 → 第二列是 。合起来：
于是从初始向量 出发，右乘一次 前进一项，乘 次就前进 项。 把 次乘法凑成一个整体，就是 ——「跨过一大段递推」被压缩成「求一个矩阵的高次幂」。可以验证 的左上角恰是 。
本质 · 递推 ↔ 矩阵，一次乘法 = 一步递推
线性递推与矩阵乘法是同一件事的两种写法：把「本步要用到的历史项」摆成状态向量，转移的每个系数就落进矩阵的一列；右乘一次 = 递推推进一步，右乘 = 一口气推进 步。 原本 的逐项递推，就此转化为「求矩阵幂 」这个能被快速幂加速的问题。
快速幂：把 Mⁿ 的 n 次拆成 log n 层
剩下的问题是怎样快速求 。若老老实实 连乘，还是 次矩阵乘法，白忙一场。快速幂（二进制倍增）的思路：把指数 写成二进制，例如 ，于是
而 这串倍增幂，每个都是前一个平方得来（），只需 次平方就能全部算出。 再按 的二进制里哪几位是 1，把对应的倍增幂累乘进结果即可。
总代价： 次平方 + 至多 次累乘，共约 次矩阵乘法；每次矩阵乘法是 （ 为矩阵阶数）。 合起来 。斐波那契 、 时也不过约 层倍增、百余次 乘法，瞬间出解。把这套骨架写成中文伪代码：
# 矩阵快速幂：求 base^n（base 是 k×k 矩阵）
res = 单位阵 I # 任何矩阵乘 I 不变，作累乘的起点
while n > 0:
if n & 1 == 1: # 当前二进制最低位为 1
res = res × base # ★把这一位对应的倍增幂累乘进结果
base = base × base # 平方 → 得到下一个倍增幂 M^(2^k)
n >>= 1 # 右移一位，看下一位
return res # 循环 ⌊log₂ n⌋+1 次，每次一到两回矩阵乘法看 Mⁿ 在倍增里累乘出来
深化 · 非标准递推：怎么构造转移矩阵
斐波那契的 太经典，容易让人以为矩阵是「背下来」的。真正的功夫在于面对一个陌生递推，自己把矩阵拼出来。 来看一个非标准递推（正是例题 P1939 的模型）：——新项跳过了 ，只由 和 决定。
第一步，定状态向量。新项用到往前三项（最远到 ），所以状态要装下最近三项：，是 3 维向量，矩阵就是 。 推进一步后，新状态应当是 。
第二步，逐行确定系数——盯住新状态的每个分量由旧分量怎样组合：
读这三行——第一行（算 ）：，系数就是递推式里的系数，写成 ——这是唯一「真正做递推加法」的一行。第二行（算 ）：新状态里的 不过是旧状态里现成的 ，原样搬过来 → 。第三行（算 ）：同理搬旧的 → 。 后两行统称位移行，作用只是把向量整体「下移一格」，好腾出位置给新项——它们对任何这类递推都长得差不多，套路固定。
构造要诀 · 一行递推系数，其余全是位移
面对形如 的常系数线性递推：状态向量取最近的若干项（维数 = 递推用到的最远回溯步数）； 转移矩阵第一行直接填递推系数 ，其余每一行是把旧分量原样下移的位移行（一个 1，其余 0）。 矩阵一旦搭好，求第 项就是算 乘初始向量，复杂度 ——递推越是「奇形怪状」，这套矩阵化越显威力。
看转移矩阵怎么从递推拼出来
a[x] 由递推给出 = a[x-1] + a[x-3]，故本行系数就是递推系数。
例题
#include <iostream>
using namespace std;
typedef long long ll;
const ll MOD = 1000000007;
// 2×2 矩阵，封装乘法（每步取模，防溢出）。斐波那契转移 M = {{1,1},{1,0}}。
struct Mat
{
ll a[2][2];
};
Mat mul(const Mat &x, const Mat &y) // 矩阵乘法：C[i][j] = Σ x[i][k]·y[k][j]
{
Mat c;
for (int i = 0; i < 2; i++)
for (int j = 0; j < 2; j++)
{
ll s = 0;
for (int k = 0; k < 2; k++)
{
s = (s + x.a[i][k] * y.a[k][j]) % MOD; // ★显式 long long + 逐步取模
}
c.a[i][j] = s;
}
return c;
}
Mat power(Mat base, ll n) // 矩阵快速幂：base^n
{
Mat res = {{{1, 0}, {0, 1}}}; // 单位阵起步
while (n > 0)
{
if (n & 1) // 当前二进制位为 1 → 累乘
{
res = mul(res, base);
}
base = mul(base, base); // 平方倍增
n >>= 1;
}
return res;
}
int main()
{
ll n;
cin >> n;
Mat M = {{{1, 1}, {1, 0}}};
Mat r = power(M, n - 1); // F(n) = (M^{n-1})[0][0]，F(1)=F(2)=1
cout << r.a[0][0] << endl;
return 0;
}
// TAG: 矩阵DP 矩阵快速幂 斐波那契 取模#include <iostream>
using namespace std;
typedef long long ll;
const ll MOD = 1000000007;
int n; // 矩阵阶数（n×n）
struct Mat
{
ll a[105][105];
};
Mat mul(const Mat &x, const Mat &y) // 通用 n×n 矩阵乘法，O(n³)
{
Mat c;
for (int i = 1; i <= n; i++)
for (int j = 1; j <= n; j++)
{
ll s = 0;
for (int k = 1; k <= n; k++)
{
s = (s + x.a[i][k] * y.a[k][j]) % MOD;
}
c.a[i][j] = s;
}
return c;
}
int main()
{
ll k;
cin >> n >> k;
Mat A, res;
for (int i = 1; i <= n; i++) // 读入待幂的矩阵
for (int j = 1; j <= n; j++)
{
cin >> A.a[i][j];
}
for (int i = 1; i <= n; i++) // res 初始化为单位阵
for (int j = 1; j <= n; j++)
{
res.a[i][j] = (i == j) ? 1 : 0;
}
while (k > 0) // 快速幂：A^k
{
if (k & 1)
{
res = mul(res, A);
}
A = mul(A, A);
k >>= 1;
}
for (int i = 1; i <= n; i++) // 输出结果矩阵
{
for (int j = 1; j <= n; j++)
{
cout << res.a[i][j] << " \n"[j == n];
}
}
return 0;
}
// TAG: 矩阵DP 矩阵快速幂 模板#include <iostream>
#include <cstring>
using namespace std;
typedef long long ll;
const ll MOD = 1000000007; // 本题模 1e9+7
// a[x] = a[x-1] + a[x-3]，a[1]=a[2]=a[3]=1。
// 状态向量 (a[x-1], a[x-2], a[x-3]) → 转移矩阵 {{1,0,1},{1,0,0},{0,1,0}}。
struct Mat
{
ll a[3][3];
};
Mat mul(const Mat &x, const Mat &y)
{
Mat c;
memset(c.a, 0, sizeof(c.a));
for (int i = 0; i < 3; i++)
for (int k = 0; k < 3; k++)
{
if (x.a[i][k] == 0) continue; // 稀疏小优化，可省
for (int j = 0; j < 3; j++)
{
c.a[i][j] = (c.a[i][j] + x.a[i][k] * y.a[k][j]) % MOD;
}
}
return c;
}
Mat power(Mat base, ll n)
{
Mat res;
memset(res.a, 0, sizeof(res.a));
for (int i = 0; i < 3; i++) res.a[i][i] = 1; // 单位阵
while (n > 0)
{
if (n & 1) res = mul(res, base);
base = mul(base, base);
n >>= 1;
}
return res;
}
int main()
{
int T;
cin >> T;
Mat M = {{{1, 0, 1}, {1, 0, 0}, {0, 1, 0}}};
while (T--)
{
ll n;
cin >> n;
if (n <= 3) // 前三项直接答 1
{
cout << 1 << endl;
continue;
}
Mat r = power(M, n - 3); // 从 (a3,a2,a1) 推到第 n 项
// 结果向量第 0 分量 = a[n] = r 作用在初始向量 (1,1,1) 上的首行之和
ll ans = (r.a[0][0] + r.a[0][1] + r.a[0][2]) % MOD;
cout << ans << endl;
}
return 0;
}
// TAG: 矩阵DP 矩阵加速 自构转移矩阵 取模