「学习笔记」矩阵
本文部分内容来自 \(\texttt{OI-Wiki}\)。
定义
对于矩阵 \(A\)，主对角线是指 \(A_{i,i}\) 的元素。
\[A =
\begin{bmatrix}
a_{1, 1} & a_{1, 2} & a_{1, 3} & \cdots & a_{1, m}\\
a_{2, 1} & a_{2, 2} & a_{2, 3} & \cdots & a_{2, m}\\
a_{3, 1} & a_{3, 2} & a_{3, 3} & \cdots & a_{3, m}\\
\cdots & \cdots & \cdots & \cdots & \cdots\\
a_{n, 1} & a_{n, 2} & a_{n, 3} & \cdots & a_{n, m}\\
\end{bmatrix}
\]
一般用 \(I\) 来表示单位矩阵，就是主对角线上为 \(1\)，其余位置为 \(0\)，如下：
\[I =
\begin{bmatrix}
1 & 0 & 0 & 0 & 0\\
0 & 1 & 0 & 0 & 0\\
0 & 0 & 1 & 0 & 0\\
0 & 0 & 0 & 1 & 0\\
0 & 0 & 0 & 0 & 1\\
\end{bmatrix}
\]
在代码中，\(n\) 行 \(m\) 列的二维数组就是一个 \(n\) 行 \(m\) 列的矩阵。
运算
矩阵加法
对于两个等行等列的矩阵，把对应位置上的数加起来就行了，如下：
\[\begin{bmatrix}
1 & 2 & 3\\
4 & 5 & 6\\
\end{bmatrix}
+
\begin{bmatrix}
2 & 3 & 3\\
2 & 3 & 3\\
\end{bmatrix}
=
\begin{bmatrix}
3 & 5 & 6\\
6 & 8 & 9\\
\end{bmatrix}
\]
矩阵减法
对于两个等行等列的矩阵，对应位置上的数相减，如下：
\[\begin{bmatrix}
1 & 2 & 3\\
4 & 5 & 6\\
\end{bmatrix}
-
\begin{bmatrix}
2 & 3 & 3\\
2 & 3 & 3\\
\end{bmatrix}
=
\begin{bmatrix}
-1 & -1 & 0\\
2 & 2 & 3\\
\end{bmatrix}
\]
数乘矩阵
类似于高中数学中的数乘向量，即将矩阵中的每一个元素都乘上这个数，如下：
\[2 \cdot
\begin{bmatrix}
1 & 2 & 3\\
4 & 5 & 6\\
\end{bmatrix}
=
\begin{bmatrix}
2 & 4 & 6\\
8 & 10 & 12\\
\end{bmatrix}
\]
矩阵乘法
矩阵 \(A\) 能与矩阵 \(B\) 做乘法，当且仅当 \(A\) 是 \(n\) 行 \(m\) 列，\(B\) 是 \(m\) 行 \(k\) 列。
（只有当第一个矩阵的列数等于第二个矩阵的行数是才能进行矩阵乘法，结果矩阵的行数是第一个矩阵的行数，列是第二个矩阵的列数）
设矩阵 \(C\) 为矩阵 \(A\) 与 \(B\) 的乘积，其中矩阵 \(C\) 中的第 \(i\) 行第 \(j\) 列元素可以表示为：
\[C_{i,j} = \sum_{k=1}^m A_{i,k}B_{k,j}
\]
如果看不懂，看一下这个例子：
\[\begin{bmatrix}
1 & 2\\
3 & 4\\
\end{bmatrix}
\times
\begin{bmatrix}
1 & 2 & 3\\
2 & 3 & 3\\
\end{bmatrix}
=
\begin{bmatrix}
5 & 8 & 9\\
11 & 18 & 21\\
\end{bmatrix}
\]
\[C_{1, 1} = A_{1, 1} \times B_{1, 1} + A_{1, 2} \times B_{2, 1}\\
C_{1, 2} = A_{1, 1} \times B_{2, 1} + A_{1, 2} \times B_{2, 2}\\
C_{1, 3} = A_{1, 1} \times B_{3, 1} + A_{1, 2} \times B_{2, 3}\\
\cdots
\]
矩阵乘法在 OI 中的应用是很多的，让我们先用代码来实现我们的矩阵乘法。
struct matrix {
int n, m;
int z[10][10];
matrix() { // 构造函数
n = m = 0;
memset(z, 0, sizeof(z));
}
};
matrix m1, m2, m3;
matrix operator * (const matrix &m1, const matrix &m2) { // 重载运算符
matrix m;
m.n = m1.n, m.m = m2.m;
for (int i = 1; i <= m.n; ++ i) {
for (int k = 1; k <= m1.m; ++ k) {
for (int j = 1; j <= m.m; ++ j) {
m.z[i][j] += m1.z[i][k] * m2.z[k][j];
}
}
}
return m;
}
m3 = m1 * m2;
矩阵乘法没有交换律，但是有结合律。
题目
P3390 【模板】矩阵快速幂
真·模板题
/*
The code was written by yifan, and yifan is neutral!!!
*/
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
template<typename T>
inline T read() {
T x = 0;
bool fg = 0;
char ch = getchar();
while (ch < '0' || ch > '9') {
fg |= (ch == '-');
ch = getchar();
}
while (ch >= '0' && ch <= '9') {
x = (x << 3) + (x << 1) + (ch ^ 48);
ch = getchar();
}
return fg ? ~x + 1 : x;
}
const int mod = 1e9 + 7;
ll k, n;
struct mat {
int h, l;
ll a[110][110];
mat() {
h = l = 0;
memset(a, 0, sizeof a);
}
} g, x;
mat operator * (const mat& x, const mat& y) {
mat res;
int h = x.h, l = y.l;
for (int i = 1; i <= h; ++ i) {
for (int k = 1; k <= l; ++ k) {
for (int j = 1; j <= l; ++ j) {
res.a[i][j] = (res.a[i][j] + (x.a[i][k] * y.a[k][j] % mod)) % mod;
}
}
}
res.h = h, res.l = l;
return res;
}
mat operator ^ (mat bas, ll x) {
mat res;
for (int i = 1; i <= n; ++ i) {
res.a[i][i] = 1;
}
res.h = bas.h, res.l = bas.l;
while (x) {
if (x & 1) {
res = res * bas;
}
bas = bas * bas;
x >>= 1;
}
return res;
}
int main() {
n = read<ll>(), k = read<ll>();
for (int i = 1; i <= n; ++ i) {
for (int j = 1; j <= n; ++ j) {
g.a[i][j] = read<ll>();
}
}
g.h = g.l = n;
g = g ^ k;
for (int i = 1; i <= n; ++ i) {
for (int j = 1; j <= n; ++ j) {
cout << g.a[i][j] << ' ';
}
putchar('\n');
}
return 0;
}
P1939 【模板】矩阵加速（数列）
一道很好的入门题。
构造矩阵
\[\begin{bmatrix}
a_{x - 1} & a_{x - 2} & a_{x - 3}
\end{bmatrix}
\]
要将该矩阵转化为
\[\begin{bmatrix}
a_x & a_{x - 1} & a_{x - 2}\\
\end{bmatrix}
\]
我们可以考虑构造一个矩阵，使得能够实现该转化。
由公式可知
\[a_{x} = 1 \cdot a_{x - 1} + 0 \cdot a_{x -2} + 1 \cdot a_{x - 3}\\
a_{x - 1} = 1 \cdot a_{x - 1} + 0 \cdot a_{x -2} + 0 \cdot a_{x - 3}\\
a_{x - 1} = 0 \cdot a_{x - 1} + 1 \cdot a_{x -2} + 0 \cdot a_{x - 3}\\
\]
我们能够发现
\[\begin{bmatrix}
a_{x - 1} & a_{x - 2} & a_{x - 3}
\end{bmatrix}
\times
\begin{bmatrix}
1 & 1 & 0\\
0 & 0 & 1\\
1 & 0 & 0\\
\end{bmatrix}
\Rightarrow
\begin{bmatrix}
a_x & a_{x - 1} & a_{x - 2}\\
\end{bmatrix}
\]
使用矩阵快速幂来优化加速。
/*
The code was written by yifan, and yifan is neutral!!!
*/
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
template<typename T>
inline T read() {
T x = 0;
bool fg = 0;
char ch = getchar();
while (ch < '0' || ch > '9') {
fg |= (ch == '-');
ch = getchar();
}
while (ch >= '0' && ch <= '9') {
x = (x << 3) + (x << 1) + (ch ^ 48);
ch = getchar();
}
return fg ? ~x + 1 : x;
}
const int mod = 1e9 + 7;
int T, n;
struct mat {
int h, l;
ll a[5][5];
mat() {
h = l = 0;
memset(a, 0, sizeof a);
}
void build() {
memset(a, 0, sizeof a);
}
} g, x;
mat operator * (const mat& x, const mat& y) {
mat res;
int h = x.h, l = y.l;
for (int i = 1; i <= h; ++ i) {
for (int k = 1; k <= l; ++ k) {
for (int j = 1; j <= l; ++ j) {
res.a[i][j] = (res.a[i][j] + (x.a[i][k] * y.a[k][j] % mod)) % mod;
}
}
}
res.h = h, res.l = l;
return res;
}
mat operator ^ (mat bas, ll x) {
mat res;
for (int i = 1; i <= 3; ++ i) {
res.a[i][i] = 1;
}
res.h = bas.h, res.l = bas.l;
while (x) {
if (x & 1) {
res = res * bas;
}
bas = bas * bas;
x >>= 1;
}
return res;
}
int main() {
T = read<int>();
while (T --) {
n = read<int>();
if (n <= 3) {
puts("1");
continue;
}
g.build(), x.build();
g.h = 1, g.l = 3;
g.a[1][1] = g.a[1][2] = g.a[1][3] = 1;
x.h = x.l = 3;
x.a[1][1] = x.a[1][2] = x.a[2][3] = x.a[3][1] = 1;
x = x ^ (n - 3);
g = g * x;
cout << g.a[1][1] % mod << '\n';
}
return 0;
}
P4910 帕秋莉的手环
这是一道矩阵优化 DP 的很好的入门题。
我们先考虑朴素的 DP。
状态：\(dp(i, 0/1)\) ：第 \(i\) 个点，是 \((1)\) 不是 \((0)\) 金属性。
转移：\(dp(i, 0) = dp(i - 1, 1), dp(i, 1) = dp(i - 1, 0) + dp(i - 1, 1)\)
然而，\(n \leq 10^{18}\)，先不说时间，空间也已经爆掉了。
我们往优化的方面想，每个状态的转移都只与上一状态有关，且都是加法，可以考虑矩阵乘法。
设矩阵
\[\begin{bmatrix}
dp(i, 1) & dp(i, 0)
\end{bmatrix}
\]
要将它转移为
\[\begin{bmatrix}
dp(i + 1, 1) & dp(i + 1, 0)
\end{bmatrix}
\]
矩阵转移式
\[\begin{bmatrix}
dp(i, 1) & dp(i, 0)\\
\end{bmatrix}
\times
\begin{bmatrix}
1 & 1\\
1 & 0\\
\end{bmatrix}
\Rightarrow
\begin{bmatrix}
dp(i + 1, 1) & dp(i + 1, 0)
\end{bmatrix}
\]
/*
The code was written by yifan, and yifan is neutral!!!
*/
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
template<typename T>
inline T read() {
T x = 0;
bool fg = 0;
char ch = getchar();
while (ch < '0' || ch > '9') {
fg |= (ch == '-');
ch = getchar();
}
while (ch >= '0' && ch <= '9') {
x = (x << 3) + (x << 1) + (ch ^ 48);
ch = getchar();
}
return fg ? ~x + 1 : x;
}
const int N = 1e6 + 5;
const int mod = 1e9 + 7;
struct mat {
int h, c;
ll z[3][3];
mat() {
h = 0, c = 0;
memset(z, 0, sizeof z);
}
void init() {
h = 0, c = 0;
memset(z, 0, sizeof z);
}
} g, f;
mat operator * (const mat& a, const mat& b) {
mat res;
int h = a.h, c = b.c;
for (int i = 1; i <= h; ++ i) {
for (int k = 1; k <= c; ++ k) {
for (int j = 1; j <= c; ++ j) {
res.z[i][j] = (res.z[i][j] + (a.z[i][k] * b.z[k][j]) % mod) % mod;
}
}
}
res.h = h, res.c = c;
return res;
}
mat operator ^ (mat x, ll y) {
mat ans;
ans.h = ans.c = 2;
for (int i = 1; i <= 2; ++ i) {
ans.z[i][i] = 1;
}
while (y) {
if (y & 1) {
ans = ans * x;
}
y >>= 1;
x = x * x;
}
return ans;
}
void solve() {
ll n = read<ll>(), ans = 0;
g.init(), f.init();
g.h = 1, g.c = 2;
f.h = f.c = 2;
g.z[1][1] = 1, g.z[1][2] = 0;
f.z[1][1] = f.z[1][2] = f.z[2][1] = 1;
f = f ^ (n - 1);
g = g * f;
ans = (ans + g.z[1][1] + g.z[1][2]) % mod;
g.init(), f.init();
g.h = 1, g.c = 2;
f.h = f.c = 2;
g.z[1][1] = 0, g.z[1][2] = 1;
f.z[1][1] = f.z[1][2] = f.z[2][1] = 1;
f = f ^ (n - 1);
g = g * f;
ans = (ans + g.z[1][1]) % mod;
cout << ans << '\n';
}
int main() {
int T = read<int>();
while (T --) {
solve();
}
return 0;
}
朝气蓬勃 后生可畏
