【精选】矩阵加速

合集 - 算法竞赛(13)

1.快速幂2024-06-07

2.【精选】矩阵加速2024-06-07

3.【二分答案】P2390 地标访问2024-06-094.动态规划初步2024-06-095.最短路算法之：Dijkstra 算法2024-06-096.最短路算法之：floyd 算法2024-06-097.最短路算法之：SPFA 算法2024-06-098.【习题】区间型动态规划2024-06-099.SPFA 判负环2024-06-0910.负环的习题和应用2024-06-0911.详解二分查找2024-06-0912.拓扑排序2024-06-0913.CF Round 920 (Div. 3) 笔记2024-06-09

收起

大家好，我是Weekoder！

今天要讲的内容是矩阵加速！

这时候就有人说了：

$\mathtt{\text{Weekoder 这么蒻，怎么会矩阵啊。还给我们讲，真是十恶不赦！}}$

不不不，容我解释。在经过我的研究后，我发现基本的矩阵运算和矩阵加速都并没有那么难。只要继续往下看，相信你也能学会！

注意：以下内容的学习难度将会用颜色表示，与洛谷题目难度顺序一致，即 $\mathtt{\text{红}}<\mathtt{\text{橙}}<\mathtt{\text{黄}}<\mathtt{\text{绿}}$。（并不对标洛谷题目难度，只作为学习难易度参考）

$$\mathtt{\text{Part 1}}\mathtt{\text{ Definition}}$$

$$\mathtt{\text{定义}}$$

矩阵和二维数组很像，是由 $m\times n$ 个数排列成 $m$ 行 $n$ 列的一张表，由于排列出来的表是一个矩形，故称其为矩阵。矩阵长这个样子：

$$\left(\begin{array}{ccccc}{a}_{11}& a_{12}& a_{13}& \cdots & a_{1n}\\ a_{21}& a_{22}& a_{23}& \cdots & a_{2n}\\ a_{31}& a_{32}& a_{33}& \cdots & a_{3n}\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ a_{m1}& a_{m2}& a_{m3}& \cdots & a_{mn}\end{array}\right)$$

可以看到，矩阵中的每个元素都有着对应的行和列，我们把一个矩阵记作 $A$，第 $i$ 行 $j$ 列的元素即为 $a_{ij}$。更形式化的，写作：

$$A=(a_{ij})\in {\mathbb{F}}^{\mathbb{m}\times \mathbb{n}}$$

其中 $\mathbb{F}$ 为数域，一般取为实数域 $\mathbb{R}$ 或复数域 $\mathbb{C}$。（看不懂没事，蒟蒻自行走开QWQ）

$$\mathtt{\text{Part 2}}\mathtt{\text{ Special matrices}}$$

$$\mathtt{\text{特殊矩阵}}$$

$\mathtt{\text{1.零矩阵}}$

元素全部为 $0$ 的矩阵称为零矩阵。像这样：

$$\left(\begin{array}{cccc}0& 0& \cdots & 0\\ 0& 0& \cdots & 0\\ ⋮& ⋮& \ddots & ⋮\\ 0& 0& \cdots & 0\end{array}\right)$$

零矩阵记作 $0_{m\times n}$，就是在 $0$ 下面加上矩阵的大小 $m\times n$。你可以把零矩阵看做数字 $0$，任何数乘以 $0$ 都得 $0$。

$\mathtt{\text{2.对角矩阵}}$

只有主对角线上的元素有值，其余元素为 $0$ 的矩阵称为对角矩阵。

注：主对角线为矩阵中从左上角到右下角的一条对角线。

$$\left(\begin{array}{cccc}{a}_1& 0& \cdots & 0\\ 0& a_2& \cdots & 0\\ ⋮& ⋮& \ddots & ⋮\\ 0& 0& \cdots & a_n\end{array}\right)$$

对角矩阵根据主对角线上的值，记作 $\text{diag(}{\text{a}}_1,{\text{a}}_2,\dots ,{\text{a}}_n\text{)}$。

$\mathtt{\text{3.单位矩阵}}$

主对角线上的元素均为 $1$，其余元素为 $0$ 的矩阵称为单位矩阵。

$$\left(\begin{array}{cccc}1& 0& \cdots & 0\\ 0& 1& \cdots & 0\\ ⋮& ⋮& \ddots & ⋮\\ 0& 0& \cdots & 1\end{array}\right)$$

单位矩阵记作 $I$。

记得分数中的概念分数单位吗？矩阵单位和分数单位的“地位”差不多，代表的都是最基础的，最小的独立个体。你可以把单位矩阵看做数字 $1$，任何数乘以 $1$ 都等于它本身。

最基础的，常见的特殊矩阵就是这些了。当然，还有很多的特殊矩阵，不过我们暂时用不到。

$$\mathtt{\text{Part 3}}\mathtt{\text{ Matrix operations}}$$

$$\mathtt{\text{矩阵运算}}$$

$\mathtt{\text{1.相等}}$

若对于矩阵 $A,B$，所有的 $i,j$ 都有 $a_{ij}=b_{ij}$ 且矩阵的行和列相等，则称矩阵 $A,B$ 相等。

其实就是两个矩阵长得一模一样。

$$\left(\begin{array}{ccccc}{a}_{11}& a_{12}& a_{13}& \cdots & a_{1n}\\ a_{21}& a_{22}& a_{23}& \cdots & a_{2n}\\ a_{31}& a_{32}& a_{33}& \cdots & a_{3n}\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ a_{m1}& a_{m2}& a_{m3}& \cdots & a_{mn}\end{array}\right)=\left(\begin{array}{ccccc}{b}_{11}& b_{12}& b_{13}& \cdots & b_{1n}\\ b_{21}& b_{22}& b_{23}& \cdots & b_{2n}\\ b_{31}& b_{32}& b_{33}& \cdots & b_{3n}\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ b_{m1}& b_{m2}& b_{m3}& \cdots & b_{mn}\end{array}\right)$$

$\mathtt{\text{2.矩阵加（减）法}}$

若要求 $A,B$ 两个矩阵之和，即 $C=A+B$，则对于任意 $i,j$，满足 $c_{ij}=a_{ij}+b_{ij}$。要求矩阵行列相等。

总结一句话：对应位置相加。

$$\begin{array}{rl}& \left(\begin{array}{ccccc}{a}_{11}& a_{12}& a_{13}& \cdots & a_{1n}\\ a_{21}& a_{22}& a_{23}& \cdots & a_{2n}\\ a_{31}& a_{32}& a_{33}& \cdots & a_{3n}\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ a_{m1}& a_{m2}& a_{m3}& \cdots & a_{mn}\end{array}\right)+\left(\begin{array}{ccccc}{b}_{11}& b_{12}& b_{13}& \cdots & b_{1n}\\ b_{21}& b_{22}& b_{23}& \cdots & b_{2n}\\ b_{31}& b_{32}& b_{33}& \cdots & b_{3n}\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ b_{m1}& b_{m2}& b_{m3}& \cdots & b_{mn}\end{array}\right)\\ =& \left(\begin{array}{ccccc}{a}_{11}+b_{11}& a_{12}+b_{12}& a_{13}+b_{13}& \cdots & a_{1n}+b_{1n}\\ a_{21}+b_{21}& a_{22}+b_{22}& a_{23}+b_{23}& \cdots & a_{2n}+b_{2n}\\ a_{31}+b_{31}& a_{32}+b_{32}& a_{33}+b_{33}& \cdots & a_{3n}+b_{3n}\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ a_{m1}+b_{m1}& a_{m2}+b_{m2}& a_{m3}+b_{m3}& \cdots & a_{mn}+b_{mn}\end{array}\right)\end{array}$$

矩阵加法满足交换律和结合律：

$$A+B=B+A$$

$$(A+B)+C=A+(B+C)$$

减法同理，对应位置相减。

$$\begin{array}{rl}& \left(\begin{array}{ccccc}{a}_{11}& a_{12}& a_{13}& \cdots & a_{1n}\\ a_{21}& a_{22}& a_{23}& \cdots & a_{2n}\\ a_{31}& a_{32}& a_{33}& \cdots & a_{3n}\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ a_{m1}& a_{m2}& a_{m3}& \cdots & a_{mn}\end{array}\right)-\left(\begin{array}{ccccc}{b}_{11}& b_{12}& b_{13}& \cdots & b_{1n}\\ b_{21}& b_{22}& b_{23}& \cdots & b_{2n}\\ b_{31}& b_{32}& b_{33}& \cdots & b_{3n}\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ b_{m1}& b_{m2}& b_{m3}& \cdots & b_{mn}\end{array}\right)\\ =& \left(\begin{array}{ccccc}{a}_{11}-b_{11}& a_{12}-b_{12}& a_{13}-b_{13}& \cdots & a_{1n}-b_{1n}\\ a_{21}-b_{21}& a_{22}-b_{22}& a_{23}-b_{23}& \cdots & a_{2n}-b_{2n}\\ a_{31}-b_{31}& a_{32}-b_{32}& a_{33}-b_{33}& \cdots & a_{3n}-b_{3n}\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ a_{m1}-b_{m1}& a_{m2}-b_{m2}& a_{m3}-b_{m3}& \cdots & a_{mn}-b_{mn}\end{array}\right)\end{array}$$

$\mathtt{\text{3.矩阵数乘}}$

数 $\lambda$（一个数字） 乘以矩阵 $A$，记作 $\lambda A$，即为矩阵数乘运算。若有 $B=\lambda A$，则对于任意 $i,j$ 都满足 $b_{ij}=\lambda a_{ij}$。

还是一句话：对应位置相乘。

$$\lambda \begin{array}{rl}& \left(\begin{array}{ccccc}{a}_{11}& a_{12}& a_{13}& \cdots & a_{1n}\\ a_{21}& a_{22}& a_{23}& \cdots & a_{2n}\\ a_{31}& a_{32}& a_{33}& \cdots & a_{3n}\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ a_{m1}& a_{m2}& a_{m3}& \cdots & a_{mn}\end{array}\right)=\left(\begin{array}{ccccc}\lambda a_{11}& \lambda a_{12}& \lambda a_{13}& \cdots & \lambda a_{1n}\\ \lambda a_{21}& \lambda a_{22}& \lambda a_{23}& \cdots & \lambda a_{2n}\\ \lambda a_{31}& \lambda a_{32}& \lambda a_{33}& \cdots & \lambda a_{3n}\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ \lambda a_{m1}& \lambda a_{m2}& \lambda a_{m3}& \cdots & \lambda a_{mn}\end{array}\right)\end{array}$$

$$\mathtt{\text{Part 3.5}}\mathtt{\text{ Matrix multiplication}}$$

$$\mathtt{\text{矩阵乘法}}$$

虽然矩阵乘法也属于矩阵运算，但难度比前面的都高，而且是今天的重点内容，所以单独放出来讲，故记为 $\mathtt{\text{Part 3.5}}$。（话说你们没有发现难度变成黄了吗）

上例题！（虽然难度是橙）

先看矩阵乘法的定义：若有 $n$ 行 $m$ 列矩阵 $A$ 和 $m$ 行 $k$ 列的矩阵 $B$（$A$ 的行与 $B$ 的列相等），则 $n$ 行 $k$ 列的矩阵 $C=A\times B$ 满足

$$c_{ij}=\sum _{l=1}^m{a}_{il}\times b_{lj}$$

只要枚举 $i,j$（范围是 $n,k$），并套用公式就能用 $O(n^3)$ 的时间复杂度解决这个问题。

我知道，这看起来根本不是新手蒟蒻能看懂的。那我就用人话来讲讲矩阵乘法。

矩阵乘法并不是一个一个乘，而是行对应列乘。怎么个乘法呢？我们来看看下面两个矩阵相乘的例子。

$$\left(\begin{array}{ccc}5& 2& 3\\ 7& 9& 4\end{array}\right)\times \left(\begin{array}{cccc}2& 6& 8& 1\\ 0& 9& 1& 3\\ 2& 4& 4& 1\end{array}\right)$$

第一个矩阵为 $A$，第二个矩阵为 $B$。

我们先取出 $A$ 的第一行。像这样：

$$\left(\begin{array}{ccc}5& 2& 3\end{array}\right)$$

再取出 $B$ 的第一列。像这样：

$$\left(\begin{array}{c}2\\ 0\\ 2\end{array}\right)$$

不对，你给我转过来。

$$\left(\begin{array}{ccc}2& 0& 2\end{array}\right)$$

现在终于可以相乘了。逐位相乘得出结果：

$$\left(\begin{array}{ccc}5& 2& 3\end{array}\right)\times \left(\begin{array}{ccc}2& 0& 2\end{array}\right)=\left(\begin{array}{ccc}5\times 2& 2\times 0& 3\times 2\end{array}\right)=\left(\begin{array}{ccc}10& 0& 6\end{array}\right)$$

得出了结果 $\left(\begin{array}{ccc}10& 0& 6\end{array}\right)$。再将每一位相加：

$$\left(\begin{array}{ccc}10& 0& 6\end{array}\right)\to 10+0+6=16$$

还记得我们之前是怎么取的吗？我们取了 $A$ 的第一行和 $B$ 的第一列（注意加粗部分），所以答案就存储在 $C$ 的第一行第一列。还没搞懂？更通用一点：我们取了 $A$ 的第 $x$ 行和 $B$ 的第 $y$ 列（注意加粗部分），所以答案就存储在 $C$ 的第 $x$ 行第 $y$ 列。也就是说，当我们想要获取矩阵 $C$ 的第 $x$ 行 $y$ 列的时候，就需要取 $A$ 的第 $x$ 行和 $B$ 的第 $y$ 列，相乘再相加。由于 $A$ 的行数与 $B$ 的列数相等，取出来的数列才可以逐位相乘（不然元素个数不一样）。而取出来的数列长度就是 $m$，所以可以用 $O(m)$ 求和，总时间复杂度 $O(nmk)=O(n^3)$。

最后，可以看看代码辅助理解。

#include <bits/stdc++.h>
using namespace std;

const int N = 105;

int n, m, k, a[N][N], b[N][N]; // 用二维数组存矩阵 A,B

int main() {
    cin >> n >> m >> k;
    for (int i = 1; i <= n; i++)
        for (int j = 1; j <= m; j++)
            cin >> a[i][j]; // 输入矩阵 A
    for (int i = 1; i <= m; i++)
        for (int j = 1; j <= k; j++)
            cin >> b[i][j]; // 输入矩阵 B
    for (int i = 1; i <= n; i++) {
        for (int j = 1; j <= k; j++) { // 枚举 C 矩阵 n 行 k 列的每个元素
  	        // 以下部分为模拟刚刚讲的矩阵乘法
            int sum = 0; // 求和，sum 即为 C_ij
            for (int l = 1; l <= m; l++)
                sum += a[i][l] * b[l][j]; // 求和，A 的行和 B 的列，建议模拟一下过程加强理解
            cout << sum << " "; // 输出 sum（C_ij）
        }
        cout << "\n"; // 记得换行！
    }
    return 0; // 完美的结束
} 

这样就能愉快地切掉这道题了。请完成这道题再继续！

矩阵乘法满足以下性质：

结合律：$(AB)C=A(BC)$

分配律：$(A+B)C=AC+BC$

矩阵乘法不满足交换律。（这是重点！）

有了矩阵乘法，我们还可以结合上面的特殊矩阵得到一些性质：

$$A\times I=A$$

$$A\times 0_{m\times n}=0_{m\times n}$$

$$\mathtt{\text{Part 4}}\mathtt{\text{ Matrix fast power}}$$

$$\mathtt{\text{矩阵封装 \& 矩阵快速幂}}$$

快到今天的主题了！上例题！

点开题目后的你 be like：

这是啥呀？

我来让题目描述“缩点水”：

给定一个 $n$ 行 $n$ 列的矩阵 $A$，求 $A^k$，即 $\underset{k\mathtt{\text{ 次}}}{\underbrace{A\times A\times A\times \cdots \times A\times A}}$。

第一思路：暴力！直接做 $k$ 次矩阵乘法，时间复杂度 $O(k{n}^3)$。看看数据范围：

$0\le k\le {10}^{12}$

考虑放弃做题。

那我们该怎么优化呢？看到需要计算 $A^k$，我突然想到了一个算法：快速幂！但是矩阵快速幂该怎么写呢？答案是：和正常的快速幂一样，矩阵也能使用快速幂，只不过快速幂中的乘法变成了矩阵乘法。但是矩阵乘法太难写，有没有什么办法能让矩阵乘法也像普通的乘法一样，只要写一个 * 乘号就行了呢？

注意：不会快速幂的话可以先简单看看我写的文章。

回到主题，有没有什么办法能只要写一个 * 乘号就能进行矩阵乘法呢？其实我们可以用结构体把矩阵封装起来，再用重载运算符就行了。关于重载运算符，可以参考这些资料。

定义一个矩阵类型的结构体可以写成这样：

struct Matrix {
	
};

我们需要在里面用一个二维数组存储矩阵。我们还可以写一个结构体初始化函数，只要定义了一个矩阵，就自动清零，免去清零的麻烦。

struct Matrix {
	int a[N][N]; // N 为矩阵大小
	Matrix() {
		memset(a, 0, sizeof a);
	}
};

最后，把矩阵乘法写进去。

struct Matrix {
    ll a[N][N];
    Matrix() {
        memset(a, 0, sizeof a);
    }
    Matrix operator*(const Matrix &x)const {
        Matrix res;
        for (int i = 1; i <= n; i++)
            for (int j = 1; j <= n; j++)
                for (int k = 1; k <= n; k++)
                    res.a[i][j] = (res.a[i][j] % MOD + a[i][k] % MOD * x.a[k][j] % MOD) % MOD;
        return res;
    }
}; 

注意，这里一定要写成 a[i][k] * x.a[k][j]，不能写成 x.a[i][k] * a[k][j]，因为矩阵乘法不满足交换律！

这样，结构体封装部分就完成了。

我们要定义两个矩阵：$a$ 和 $base$。$a$ 是输入的矩阵，$base$ 是答案矩阵，所以 $base$ 需要初始化成 $I$（单位矩阵），写一个初始化函数 $\mathrm{init}$，如下：

void init() {
    for (int i = 1; i <= n; i++) base.a[i][i] =1;
}

初始化完以后，就可以执行快速幂了，计算 $A^k$ 了，让 $base$ 乘 $A$。矩阵快速幂核心代码如下：

void expow(ll b) {
    while (b) {
        if (b & 1) base = base * a;
        a = a * a, b >>= 1;
    }  
}

有一点需要注意的就是，不能写成 base *= a 等形式，因为重载运算符定义的是 *，没有定义 *=，所以需要将 *= 展开。

最后，就可以输出 $base$ 了。展示全部代码：

#include <bits/stdc++.h>
using namespace std;

typedef long long ll;

const int N = 105, MOD = 1e9 + 7;

int n;
ll k;

struct Matrix {
    ll a[N][N];
    Matrix() {
        memset(a, 0, sizeof a);
    }
    Matrix operator*(const Matrix &x)const {
        Matrix res;
        for (int i = 1; i <= n; i++)
            for (int j = 1; j <= n; j++)
                for (int k = 1; k <= n; k++)
                    res.a[i][j] = (res.a[i][j] % MOD + a[i][k] % MOD * x.a[k][j] % MOD) % MOD;
        return res;
    }
}a, base; 

void init() {
    for (int i = 1; i <= n; i++) base.a[i][i] =1;
}

void expow(ll b) {
    while (b) {
        if (b & 1) base = base * a;
        a = a * a, b >>= 1;
    }  
}

int main() {
    cin >> n >> k;
    for (int i = 1; i <= n; i++)
        for (int j = 1; j <= n; j++)
            cin >> a.a[i][j];
    init();
    expow(k);
    for (int i = 1; i <= n; putchar('\n'), i++)
        for (int j = 1; j <= n; j++)
            cout << base.a[i][j] << " ";
    return 0;
} 

$$\mathtt{\text{Part 5}}\mathtt{\text{ Matrix acceleration}}$$

$$\mathtt{\text{矩阵加速}}$$

终于到了最后的 $\mathtt{\text{BOSS 关卡}}$ 了！你们有信心吗？加油！

点击此处进入 $\mathtt{\text{BOSS 关卡}}$ ......

点开题目 $\mathtt{\text{BOSS 关卡}}$ 后的你 be like（梅开二度）：

这和矩阵有什么关系吗？？？

我直接一个递推！

• 对于 $100\mathrm{\%}$ 的数据 $1\le T\le 100$，$1\le n\le 2\times {10}^9$。

$O(Tn)$ 这 $2\times {10}^{11}$ 的复杂度实在无法接受。

（呜呜呜我再也不学 c艹 了）

没关系，先看看思路！

因为发现当 $x\le 3$ 时答案为 $1$，所以这是最基础的情况。我们可以构造一个只有一列的矩阵：

$$\left(\begin{array}{ccc}{a}_3& a_2& a_1\end{array}\right)=\left(\begin{array}{ccc}1& 1& 1\end{array}\right)$$

显然，这三个元素都是 $1$。

那么，假设我想要得到 $a_4$，该怎么办呢？所以，我们需要进行一种运算，让上面的矩阵变化一下，像这样：

$$\left(\begin{array}{ccc}{a}_3& a_2& a_1\end{array}\right)\to \left(\begin{array}{ccc}{a}_4& a_3& a_2\end{array}\right)$$

更加通用一点：

$$\left(\begin{array}{ccc}{a}_x& a_{x-1}& a_{x-2}\end{array}\right)\to \left(\begin{array}{ccc}{a}_{x+1}& a_x& a_{x-1}\end{array}\right)$$

可以发现，矩阵中的每个元素的项数都向前推进了 $1$。那么，我们大概可以写出伪代码：

---

如果 $x\le 3$

输出 $1$

否则

执行运算 $n-3$ 次（重要！）

并输出答案矩阵 $1$ 行 $1$ 列

---

特判（对于特殊情况的判断）和输出应该没什么问题，主要是为什么运算恰好要执行 $n-3$ 次呢？稍微画个图模拟一下就好了。

还是假设要获取 $a_4$，则执行运算 $4-3=1$ 次。在执行 $1$ 次运算后，

$$\left(\begin{array}{ccc}{a}_3& a_2& a_1\end{array}\right)$$

变为

$$\left(\begin{array}{ccc}{a}_4& a_3& a_2\end{array}\right)$$

这样就刚好在第 $1$ 行 $1$ 列得到 $a_4$ 啦！

那么，说了这么久，这个神秘的运算是什么呢？当当当当~，他就是我们的——矩阵乘法！

没错，所谓的变换，其实就是乘上了一个特殊的矩阵！那么，这个矩阵长什么样呢？让我们一起来推理吧。

（此处应配上推理の小曲）

我们可以先列一个表格，表格的行代表矩阵 $\left(\begin{array}{ccc}{a}_3& a_2& a_1\end{array}\right)$ 的元素，列代表递推时与这些元素相关的元素。像这样：（表格可能在博客里渲染不出来，凑合着看吧，抱歉）

	$a_x$	$a_{x-1}$	$a_{x-2}$
$a_{x-1}$
$a_{x-2}$
$a_{x-3}$

好了，对于 $a_x$，我们该怎么填他那一列呢？我们可以观察到递推式 $a_x=a_{x-1}+a_{x-3}$，所以有：

$$a_x=a_{x-1}\times 1+a_{x-2}\times 0+a_{x-3}\times 1$$

观察系数 $1,0,1$，把这些系数填入表格中：

	$a_x$	$a_{x-1}$	$a_{x-2}$
$a_{x-1}$	$1$
$a_{x-2}$	$0$
$a_{x-3}$	$1$

后面的也以此类推：

$$a_{x-1}=a_{x-1}\times 1+a_{x-2}\times 0+a_{x-3}\times 0$$

$$a_{x-2}=a_{x-1}\times 0+a_{x-2}\times 1+a_{x-3}\times 0$$

	$a_x$	$a_{x-1}$	$a_{x-2}$
$a_{x-1}$	$1$	$1$	$0$
$a_{x-2}$	$0$	$0$	$1$
$a_{x-3}$	$1$	$0$	$0$

这样，我们就可以推出这个神秘的矩阵了：

$$\left(\begin{array}{ccc}1& 1& 0\\ 0& 0& 1\\ 1& 0& 0\end{array}\right)$$

好了，现在我们终于知道了，一次神秘操作，就是将让 $\left(\begin{array}{ccc}{a}_3& a_2& a_1\end{array}\right)$ 这个矩阵乘上$\left(\begin{array}{ccc}1& 1& 0\\ 0& 0& 1\\ 1& 0& 0\end{array}\right)$。这时候就有人问了：

一次矩阵乘法的时间复杂度还没有递推快，这根本就没有优化嘛。

等等！我们把这个式子展开：

$$\begin{array}{rl}& \left(\begin{array}{ccc}1& 1& 1\end{array}\right)\times \left(\begin{array}{ccc}1& 1& 0\\ 0& 0& 1\\ 1& 0& 0\end{array}\right)\times \left(\begin{array}{ccc}1& 1& 0\\ 0& 0& 1\\ 1& 0& 0\end{array}\right)\cdots \times \left(\begin{array}{ccc}1& 1& 0\\ 0& 0& 1\\ 1& 0& 0\end{array}\right)\\ =& \left(\begin{array}{ccc}1& 1& 1\end{array}\right)\times {\left(\begin{array}{ccc}1& 1& 0\\ 0& 0& 1\\ 1& 0& 0\end{array}\right)}^{n-3}\end{array}$$

不是吧！这居然变成了一个矩阵快速幂？！！

也就是说，我们可以用快速幂计算 ${\left(\begin{array}{ccc}1& 1& 0\\ 0& 0& 1\\ 1& 0& 0\end{array}\right)}^{n-3}$，并乘上初始矩阵 $\left(\begin{array}{ccc}1& 1& 1\end{array}\right)$。这样，我们成功地把时间复杂度从 $O(Tn)$ 优化到了 $O(T\log n)$！（矩阵快速幂是 $O(\log n)$，因为矩阵很小，矩阵乘法只计算 $9$ 次，是一个很小的常数）

下面奉上代码：（标准的矩阵加速思想）

#include <bits/stdc++.h>
using namespace std;

typedef long long ll;

const int MOD = 1e9 + 7;

int T, n;

struct Matrix {
    ll a[5][5];
    Matrix() {
        memset(a, 0, sizeof a);
    }
    Matrix operator*(const Matrix &x)const { // 矩阵乘法
        Matrix res;
        for (int i = 1; i <= 3; i++)
            for (int j = 1; j <= 3; j++)
            	for (int k = 1; k <= 3; k++)
                    res.a[i][j] = (res.a[i][j] % MOD + a[i][k] % MOD * x.a[k][j] % MOD) % MOD;
        return res;
    }
    void mems() {
    	memset(a, 0, sizeof a);
	}
}ans, base; 

void init() { // 初始化两个矩阵
	ans.mems(), base.mems(); // 记得清空！
	ans.a[1][1] = ans.a[1][2] = ans.a[1][3] = 1;
	base.a[1][1] = base.a[1][2] = base.a[2][3] = base.a[3][1] = 1;
}

void expow(int b) { // 矩阵快速幂，是在 ans 矩阵的基础上乘的
    while (b) {
        if (b & 1) ans = ans * base;
        base = base * base, b >>= 1;
    }  
}

int main() {    
    cin >> T;
    while (T --) {
        cin >> n;
        init(); // 初始化不能忘
        if (n <= 3) { // 特判
            cout << "1\n";
            continue;
        } 
        expow(n - 3); // 计算特殊矩阵的 n - 3 次方，已经乘到了 ans 里
        cout << ans.a[1][1] << "\n"; // 输出答案！芜湖！
    }
    return 0; // 快乐结束
} 

就这样，我们完成了矩阵加速递推。

再次声明矩阵快速幂（矩阵加速）时间复杂度：$O(N^3\log n)$，其中 $N$ 为矩阵的行数（列数），$n$ 为快速幂的规模 $a^n$。

小提示：关于 $base$ 矩阵的构造

就是这个 $\left(\begin{array}{ccc}1& 1& 0\\ 0& 0& 1\\ 1& 0& 0\end{array}\right)$ 矩阵。

可以这样：我们要推导出 $\left(\begin{array}{ccc}{a}_x& a_{x-1}& a_{x-2}\end{array}\right)$，那么这个矩阵从哪里来？当然是从 $\left(\begin{array}{ccc}{a}_{x-1}& a_{x-2}& a_{x-3}\end{array}\right)$ 来。所以，表格才长这样：

	$a_x$	$a_{x-1}$	$a_{x-2}$
$a_{x-1}$
$a_{x-2}$
$a_{x-3}$

那么，能不能构造一个行列数各不相同的矩阵，而不是一个 $n\times n$ 的矩阵呢？答案是不可以，因为我们要计算 $\left(\begin{array}{ccc}1& 1& 0\\ 0& 0& 1\\ 1& 0& 0\end{array}\right)$ 这种矩阵的幂，那如果行和列不相等，相乘的两个矩阵的行列也不相等，就无法进行矩阵乘法。比如这个：

$$\left(\begin{array}{ccc}1& 1& 0\\ 0& 0& 1\end{array}\right)\times \left(\begin{array}{ccc}1& 1& 0\\ 0& 0& 1\end{array}\right)$$

可以看到，左边 $2$ 行，右边 $3$ 列，显然不相等，无法进行矩阵乘法。

$$\mathtt{\text{Part 6}}\mathtt{\text{ Thank you!}}$$

$$\mathtt{\text{你}}\mathtt{\text{居}}\mathtt{\text{然}}\mathtt{\text{看}}\mathtt{\text{完}}\mathtt{\text{了}}\mathtt{\text{！}}$$

这篇文章花费了我很多时间，希望你喜欢！

对了，你学会了吗？是不是，矩阵也并没有那么难？

这应该是我的【精选】文章中的第一篇，没想到写的是矩阵方面的。

总之，很感谢你的阅读！希望你能从我这学到点东西！

再见！