学习笔记:矩阵乘法,高斯消元,线性基
矩阵乘法
引入
如果 \(C = AB\),则 \(c_{ij} = \sum\limits_{k = 1}^{n}a_{ik} \cdot b_{kj}\),即 \(A\) 的第 \(i\) 行与 \(B\) 的第 \(j\) 列的点积。
-
假设有 \(n\) 个地点,\(i\) 到 \(j\) 做飞机有 \(a_{ij}\) 种选择,坐火车有 \(b_{ij}\) 种选择。求从 \(i\) 先做飞机再坐火车到 \(j\) 的方案。
枚举中转点 \(k\)。
\(c_{ij} = \sum\limits_{k = 1}^na_{ik}\cdot b_{kj}\),写成矩乘即 \(C = AB\)。
矩阵快速幂优化线性递推
-
\(f(n) = a_1f(n - 1) + a_2f(n - 2) + a_3f(n - 3)\).
\[\begin{pmatrix} f_n\\ f_{n - 1}\\ f_{n - 2}\\ \end{pmatrix} = \begin{pmatrix} a_1 & a_2 & a_3\\ 1 & 0 & 0 \\ 0 & 1 & 0 \\ \end{pmatrix} \begin{pmatrix} f_{n - 1}\\ f_{n - 2}\\ f_{n - 3}\\ \end{pmatrix} \] -
\(f(n) = a_1f(n - 1) + a_2f(n - 2) + C\).
\[\begin{pmatrix} f_n\\ f_{n - 1}\\ C \end{pmatrix} = \begin{pmatrix} a_1 & a_2 & 1\\ 1 & 0 & 0\\ 0 & 0 & 1 \end{pmatrix} \begin{pmatrix} f_{n - 1}\\ f_{n - 2}\\ C \end{pmatrix} \]或
\[\begin{pmatrix} f_n\\ f_{n - 1}\\ 1 \end{pmatrix} = \begin{pmatrix} a_1 & a_2 & C\\ 1 & 0 & 0\\ 0 & 0 & 1 \end{pmatrix} \begin{pmatrix} f_{n - 1}\\ f_{n - 2}\\ 1 \end{pmatrix} \] -
\(f(n) = a_1f(n - 1) + a_2f(n - 2) + c_2n^2 + c_1n + c_0\).
\[\begin{pmatrix} f_n\\ f_{n - 1}\\ (n + 1)^2\\ n + 1\\ 1 \end{pmatrix} = \begin{pmatrix} a_1 & a_2 & c_2 & c_1 & c_0\\ 1 & 0 & 0 & 0 & 0 \\ 0 & 0 & 1 & 2 & 1 \\ 0 & 0 & 0 & 1 & 1 \\ 0 & 0 & 0 & 0 & 1 \\ \end{pmatrix} \begin{pmatrix} f_{n - 1}\\ f_{n - 2}\\ n^2\\ n\\ 1 \end{pmatrix} \] -
\(f(n) = a_{11}f(n - 1) + a_{12}g(n - 1), \ g(n) = a_{21}f(n - 1) + a_{22}g(n - 2)\).
\[\begin{pmatrix} f_n\\ g_n\\ \end{pmatrix} = \begin{pmatrix} a_{11} & a_{12}\\ a_{21} & a_{22}\\ \end{pmatrix} \begin{pmatrix} f_{n - 1}\\ g_{n - 1}\\ \end{pmatrix} \] -
若 \(f(n)\) 与 \(g(n)\) 是线性递推的,则 \(h(n) = f(n)\cdot g(n)\) 也能线性递推。
\[\begin{aligned} &\\ f(n) &= f(n - 1) + f(n - 2)\\ &\\ f^2(n) &= f^2(n - 1) + f^2(n - 2) + 2\cdot f(n - 1) \cdot f(n - 2)\\ &\\ f(n - 1) \cdot f(n - 2) &= \begin{cases} f(n - 1)\cdot (f(n - 1) - f(n - 3)) \\ (f(n - 2) + f(n - 3))\cdot f(n - 2) \end{cases}\\ &\\ 2\cdot f(n - 1) \cdot f(n - 2) &= f^2(n - 1) + f^2(n - 2) - f^2(n - 3)\\ &\\ f^2(n) &= 2\cdot f^2(n - 1) + 2\cdot f^2(n - 2) - f^2(n - 3) \end{aligned} \]
例题
CF1117D
题意:有多少长度为 \(n\) 的 \(01\) 串,使得任意极长全 \(0\) 字串的长度都被 \(m\) 整除?(\(n \le 10^{18}, \ m \le 100\))
定义 \(f_{i, 0/1}\) 为长度为 \(i\),结尾为 \(0/1\) 的方案数。
上下相加得到 \(f_i = f_{i - 1} + f_{i - m}\)。
构造矩阵乘法加速递推。
于是有 \(F_{n} = AF_{n - 1} = A^2F_{n - 2} = A^{n - m}F_m\)。
预处理 \(F_m\),时间复杂度 \(O(m + m^3\log n)\)。
注意特判 \(n < m\) 的情况。
NC17890
题意:\(m \times n\) 的方格,每个格子可以填白色或黑色,左右两个不同为白,左右两列不同为全黑,求合法方案数。(\(n\le10^{18}, m \le5\))
先朴素动态规划。
\(f_{i, s}\) 表示第 \(i\) 列状态为 \(s\) 时的方案,'0' 表示白,'1' 表示黑。
考虑 \(s\) 能从哪些 \(s'\) 转移来。
- \(s\ |\ s' = 2^m - 1\)。
- 若 \(s = 2^m - 1\),\(s'\ne 2^m - 1\)。
所以
又可以看成 \(F_n = A\cdot F_{n - 1} = A^{n - 1}\cdot F_1\) 的形式。
CF718C
题意:给定数列 \(\{a\}\),支持两种操作。
- 区间 \([l, r]\) 加 \(x\)。
- 求 \(\sum_{i = l}^rf(a_i)\),\(f(1) = f(2) = 1, \ f(n) = f(n - 1) + f(n - 2)\)。
写成向量的形式
根据递推式,可以写出 \(\vec F(0) = \begin{pmatrix}1\\0\end{pmatrix}\)。
于是可以用 \(\vec{F(a_i)}_{1, 1}\) 来替代原序列,将操作转化为
- 区间 \([l, r]\) 左乘 \(\begin{pmatrix} 1 & 1\\ 1 & 0\\ \end{pmatrix}^x\)。
- 查询 \(\sum_{i = l}^r\vec{F(a_i)}_{1, 1}\)。
用一个 \(2 \times 1\) 的向量表示线段树的节点,\(2 \times 2\) 的向量表示懒标记,维护区间矩阵乘法。
P6327 区间加区间 sin 和
题意:给出一个长度为 \(n\) 的整数序列 \(a_1,a_2,\ldots,a_n\),进行 \(m\) 次操作,操作分为两类。
操作 \(1\):给出 \(l,r,v\),将 \(a_l,a_{l+1},\ldots,a_r\) 分别加上 \(v\)。
操作 \(2\):给出 \(l,r\),询问 \(\sum\limits_{i=l}^{r}\sin(a_i)\)。
\(\sin(x + v)\) 和 \(\cos(x + v)\) 都是线性递推的。
线段树维护即可。
斐波*
题意:斐波那契数列 \(\{fib\}\) 满足
\(S\) 是一个可重集合 \(\left\{s_1,s_2,...,s_{|S|}\right\}\)。
\(f(S)\) 定义为
支持两种操作。
- 把 \(a_p\) 变为 \(v\)。
- 计算 \(\sum_{i=l}^{r}\sum_{j=i}^{r}f(\left\{a_i,a_{i+1},..,a_j\right\})\) 。
定义 \(g(n) = fib^2(n)\)。
证明参考上文 矩阵快速幂优化线性递推.5。
写成向量形式。
虽然 \(\vec{G(0)}\) 不存在,但仍可通过递推式求出 \(\vec{G(0)} = \begin{pmatrix}0\\1\\1\\\end{pmatrix}\)。
用向量的第一个元素表示答案。
往集合里新增一个元素 \(a\)。
其中 \(I\) 是单位矩阵。
所以
将 \((I + A^{a_i})\) 视作元素 \(v_i\)。
题目要求的 \(\sum_{i=l}^{r}\sum_{j=i}^{r}f(\left\{a_i,a_{i+1},..,a_j\right\})\) 即
即区间 \([l, r]\) 内所有子区间的元素积之和。
考虑分治。
将 \([l, r]\) 分为 \([l, mid]\) 和 \([mid + 1, r]\)。
讨论贡献来源,令 \(\{l, r\}\) 表示 \(\prod_{i = l}^r v_i\)。
- 子区间单独属于左区间或右区间,向下递归。
- 子区间跨块,贡献为 \(\sum\limits_{i = l}^{mid}\{i, mid\} \cdot \sum\limits_{i = mid + 1}^{r} \{mid + 1, i\}\).
具体实现用线段树内维护 \(4\) 个信息。
对于代表区间 \([l, r]\) 的节点。
- \(ans\) 表示 \([l, r]\) 的答案。
- \(self = \prod_{i = l}^r v_i\)
- \(left = \sum_{i = l}^r\{l, i\}\)
- \(right = \sum_{i = l}^r\{i, r\}\)
struct Node {
Matrix ans, self, l, r;
Node() { self = I; }
void operator = (Matrix x){
ans = self = l = r = x;
}
friend Node operator + (Node L, Node R) {
Node x;
x.ans = L.ans + R.ans + L.r * R.l;
x.self = L.self * R.self;
x.l = L.l + R.l * L.self;
x.r = R.r + L.r * R.self;
return x;
}
} t[N << 2];
线性方程组
Gauss-Jordan 消元法
思想很简单,即
- 选主元。
- 选主元不为零的行为主行。
- 用主行消去主行外的所有主元。
判断无解或无穷解。
如果最终增广矩阵如下。
最后一行系数矩阵为零,而增广矩阵非零,无解。
如果增广矩阵也为零呢?
可以得到
其中 \(x_2, \ x_4 \in \R\) ,有无穷多组解,\(x_2, \ x_4\) 确定,则解确定。
增广矩阵每一行第一个非零元对于的变量称为 首变量(lead variable),化简过程种跳过的变量称为 自由变量(free variable)。
因此 \(x_1, x_3, x_5\) 为首变量,\(x_2, x_4\) 为自由变量。
时间复杂度 \(O(n^3)\)。
#include<bits/stdc++.h>
using namespace std;
using db = double;
constexpr int N = 105;
db a[N][N];
int n;
bool neq(db x, db y) {
return fabs(x - y) > 1e-6;
}
int main() {
cin.tie(0)->sync_with_stdio(0);
cin >> n;
for(int i = 1; i <= n; ++ i) {
for(int j = 1; j <= n + 1; ++ j) {
cin >> a[i][j];
}
}
int cur = 1;
for(int k = 1; k <= n; ++ k) { \\ 选主元
for(int i = cur; i <= n; ++ i) {
if(neq(a[i][k], 0)) { \\ 选主行
for(int j = k; j <= n + 1; ++ j) {
swap(a[cur][j], a[i][j]);
}
break;
}
}
if(neq(a[cur][k], 0)) { \\ 消元
db x = a[cur][k];
for(int j = k; j <= n + 1; ++ j) {
a[cur][j] /= x;
}
for(int i = 1; i <= n; ++ i) {
if(i == cur) continue;
db x = a[i][k];
for(int j = k; j <= n + 1; ++ j) {
a[i][j] -= x * a[cur][j];
}
}
++ cur; \\ 主行移至下一行
}
}
for(int i = cur; i <= n; ++ i) {
if(neq(a[i][n + 1], 0)) {
cout << -1;
exit(0);
}
}
if(cur != n + 1) {
cout << 0;
}
else {
for(int i = 1; i <= n; ++ i) {
cout << "x" << i << "=" << fixed << setprecision(5) << a[i][n + 1] << '\n';
}
}
return 0;
}
异或方程组
异或可看成模 \(2\) 意义下的加法,写法与加法无异。
有几个优化的点。
bitset<N> a[N]
替代系数矩阵。- 交换用
swap(a[i], a[cur])
。 - 消元用
a[i] ^= a[cur])
。
时间复杂度 \(O(\dfrac{n^3}{w})\)。
例题
POJ1222
题意:\(5 \times 6\) 的方格,改变一盏灯的状态会改变周围所有灯的状态,给出一种开关方案使得所有灯熄灭。
令 \(x_{i, j} = 0 / 1\) 表示 \((i, j)\) 有无操作。
如果 \((i, j)\) 初状态为 \(a_{i, j}\),则达到末状态需要满足 \(a_{i, j} \oplus\)
对 \(30\) 个未知数与 \(30\) 个方程进行高斯消元。
如果处理无穷多组解的情况?将自由变量全赋值为 \(0\),第 \(i\) 行的首变量即列的增广一列的值。
CF1344F
题意:一个包含三原色 RYB
的序列混合的结果定义为:
- 如果序列开头两项颜色相同,将这两项删去。
- 如果序列开头两项颜色不同,将这两项替换为与这两种颜色不同的颜色。
- 特别地,如果序列为空,则混合的结果是白色
W
。
有一个长为 \(n\) 的颜色序列(某些位置为空),给出 \(k\) 个操作:
mix
:选择一个子序列和混合时的顺序(忽略空位置),给出其混合后的结果。RY
:选择一个子序列,将所有R
变为Y
,所有Y
变为R
,B
和空位置不变。RB
:选择一个子序列,将所有R
变为B
,所有B
变为R
,Y
和空位置不变。YB
:选择一个子序列,将所有Y
变为B
,所有B
变为Y
,R
和空位置不变。
你需要根据 mix
操作的信息确定一种可能的原序列。保证 \(n, k \leqslant 10^3\)。
将 WRYB
依次赋值为 \(0, 1, 2, 3\),可以发现 mix
操作即对取出元素做异或。
将 \((00)_2, (01)_2, (10)_2, (11)_2\) 表示为 \((b_ia_i)_2\)。
RY
,交换选定数的 \(a_i, b_i\)。RB
,\(b_i \leftarrow a_i \oplus b_i\)。YB
,\(a_i \leftarrow a_i \oplus b_i\)。
动态维护当前的 \(a_i\) 和 \(b_i\) 与初始值的关系,设 \(i\) 初始为 \((b_0a_0)_2\)。
- \(t_i = (00) \iff a_i = 0\)。
- \(t_i = (01) \iff a_i = a_{i_0}\)。
- \(t_i = (10) \iff a_i = b_{i_0}\)。
- \(t_i = (11) \iff a_i = a_{i_0} \oplus b_{i_0}\)。
\(b_i\) 同理维护一个 \(t_i'\)。
把 \(a_{i_0}\),\(b_{i_0}\) 表示为 \(2n\) 个未知数。
对于一次 mix
操作,可以对分别对当前 \(\{a_i\}\) 和\(\{b_i\}\) 异或得到关于 \(a_{i_0}\),\(b_{i_0}\)的两组方程。
高斯消元求出所有\(a_{i_0}, \ b_{i_0}\)。
线性基
引入
线性基即高中数学向量基底拓展到 \(n\) 维的情况。
线性基内元素线性无关,即不能互相表示。
线性基可以表示出任意当前线性空间内的元素。
异或线性基
对于一组基 \(\{b\}\),其内元素不能通过异或运算互相表示。
-
高斯消元求线性基
求解 \(\{a\} = {1, 3, 4, 6, 7}\) 的一组基底。
把一个 \(n\) 位二进制数看作 \(n\) 维向量。
\[\begin{bmatrix} 1 & 0 & 0\\ 1 & 1 & 0\\ 0 & 0 & 1\\ 0 & 1 & 1\\ 0 & 1 & 1\\ \end{bmatrix} \]执行高斯消元。
\[\begin{bmatrix} 1 & 0 & 0\\ 0 & 1 & 0\\ 0 & 0 & 1\\ 0 & 0 & 0\\ 0 & 0 & 0\\ \end{bmatrix} \]剩下 \(3\) 个非零元素 \(\{1, 2, 4\}\)即为 \(\{a\}\) 的一组基。
显然可用 \(1, 2, 4\) 表示出 \(3, 6, 7\) 却无法互相表示。
引理:大小为 \(sz\) 的基底能组成 \(2^{sz}\) 个不同的数。
证明:高斯消元得到的基底互不相交,任意组合互不相同。
-
在线求线性基
高斯消元是离线的,复杂且有局限性。
设 \(b_i\) 为基底中最高位为 \(i\) 的元素。
用如下方法插入新元素,其实际意义与高斯消元无异。
void insert(int x) { for(int i = 30; i >= 0; -- i) { if(x >> i & 1) { if(b[i] == 0) { b[i] = x; return; } x ^= b[i]; } } }
-
求线性基内得到的最大(小)元素
高斯消元得到的基底元素交集为空,全部异或起来即最大值。
在线插入法可以如下贪心:
int get_max() { int cur = 0; for(int i = 30; i >= 0; -- i) { cur = max(cur, cur ^ b[i]); } return cur; }
-
判断某个元素是否在线性空间内存在
低位操作不会影响高位,从高到低贪心。
bool exist(int x) { int cur = 0; for(int i = 30; i >= 0; -- i) { if((x >> i & 1) != (cur >> i & 1)) { cur ^= b[i]; } } return x == cur; }
例题
CF1101G
题意:
给定 \(\{a_i\}\),试将其划分为尽可能多的非空子段,满足每一个元素出现且仅出现在其中一个子段中,且在这些子段中任取若干段,它们包含的所有数的异或和不能为\(0\).
如果不存在方案输出 \(-1\),否则输出所有合法的划分方案中最大的划分数。
表示成前缀异或。
设 \(i_1, i_2, \dots i_k\) 为划分点。
则 \(\{s_{i_2}\oplus s_{i_1 - 1}, \dots s_{i_k}\oplus s_{i_{k - 1} - 1}\}\) 互相组合不能表示出零。
表示出零,当且仅当存在元素能被其他元素表示。
所以 \(\{s_{i_2}\oplus s_{i_1 - 1}, \dots s_{i_k}\oplus s_{i_{k - 1} - 1}\}\) 的线性基(等效于 \(\{s_i\}\) 的线性基)为极大符合条件集合。
线性基的大小唯一确定,无解的情况为 \(s_n = 0\)。