前のページ

浅学矩阵

矩阵优化基础

前情提要:开始准备明年省选了,先把 NOIp 的重要知识点梳理一遍,补一下以前没学过的东西。像这种有内容转化输出的会记录下来,私题或者流水模拟赛补题的丢在博客园里

何为矩阵

矩阵的形式化定义似乎会牵扯到 向量,然后就变成一个数学问题。

矩阵的最初目的只是为线性方程组提供一个简写的形式。

这里我们考虑肤浅地从计算机语言的角度来感性理解,矩阵并非我们经常理解的“数”,而是一个比较复杂的模型,简单记为“集合”,这个集合里封装了任意类型的值。

矩阵的加减法只需要简单地把两个矩阵对应位置上的数相加减。

$$ \boldsymbol{C} = \boldsymbol{A} \pm \boldsymbol{B} \Leftrightarrow \forall i \in [1, n], \forall j \in [1, m], \boldsymbol{C}_{i, j} = \boldsymbol{A}_{i, j} \pm \boldsymbol{B}_{i, j} $$

其中,$\boldsymbol{A, B, C}$ 都是 $n \times m$ 的矩阵。

而其在 OI 中最重要也最常见的莫过于矩阵乘法。

考虑两个矩阵的相乘:

$$ \begin{bmatrix} 2 & 1 \\ 4 & 3 \end{bmatrix} \times \begin{bmatrix} 1 & 2 \\ 1 & 0 \end{bmatrix} = \begin{bmatrix} 3 & 4 \\ 7 & 8 \end{bmatrix} $$

定义两个矩阵相乘后的“结果矩阵”中第一个矩阵第 $n$ 行与第二个矩阵第 $m$ 列 对应每个位置 的乘积之

形式化定义如下,令 $\boldsymbol{C}_{i, j}$ 表示 $\boldsymbol{A}$ 的第 $i$ 行与 $\boldsymbol{B}$ 的第 $j$ 列对应乘积的和:

$$ \boldsymbol{C}_{i, j} = \sum \boldsymbol{A}_{i, k} \times \boldsymbol{B}_{k, j} (1 \leq i \leq n, 1 \leq j \leq n, 1 \leq k \leq n) $$

参与矩阵乘法运算的第一个矩阵的列数必须等于第二个矩阵的行数,所得到的结果矩阵的行数、列数分别等于第一个矩阵的行数、第二个矩阵的列数。

关于转置等相关内容不在此赘述,后面如果线代用的多了可以专门开一个系列来记录。

注意:

  • 矩阵乘法满足结合律:$(\boldsymbol{A} \times \boldsymbol{B}) \times \boldsymbol{C} = \boldsymbol{A} \times (\boldsymbol{B} \times \boldsymbol{C})$
  • 矩阵乘法满足分配率:$(\boldsymbol{A} + \boldsymbol{B}) \times \boldsymbol{C} = \boldsymbol{A} \times \boldsymbol{C} + \boldsymbol{B} \times \boldsymbol{C}$
  • 矩阵乘法不满足交换律!

考虑一种更为特殊的情形:$\boldsymbol{F}$ 是一个 $1 \times n$ 的矩阵,$\boldsymbol{A}$ 是一个 $n \times n$ 的矩阵,则对于 $\boldsymbol{F’} = \boldsymbol{F} \times \boldsymbol{A}$ 也是 $1 \times n$ 的矩阵

把 $\boldsymbol{F, F’}$ 看做一维数组,省略其行下标 $1$,按照矩阵乘法的定义:

$$ \forall j \in [1, n], \boldsymbol{F'} = \sum_{k = 1}^{n} \boldsymbol{F}_{k} \times \boldsymbol{A}_{k, j} $$

扣住定义来理解,具体含义为,经过一次与 $\boldsymbol{A}$ 的矩阵乘法后,$\boldsymbol{F}$ 数组中的第 $k$ 个元素会以 $\boldsymbol{A}_{k, j}$ 为系数累加到 $\boldsymbol{F’}$ 的第 $j$ 个值上,这 等价于在一个向量的 $k, j$ 两个状态间发生递推

这是裸的代码,未封装:

constexpr int N = 100;

int n;
int A[N][N], B[N][N], C[N][N];

void matrixMulti() {
    memset(C, 0, sizeof(C));
    for (int i = 1; i <= n; i++) {
        for (int j = 1; j <= n; j++)
            for (int k = 1; k <= n; k++)
                C[i][j] += A[i][k] * B[k][j];
    }
}

快速幂

快速幂是一个应用了倍增思想来优化求幂时结果太大而降速的特殊乘法,核心思想就是:

$$ a^i = (a^{\frac{i}{2}})^2 $$

快速幂的实质是 不断让指数变小,度数变大从而减少运算的次数

实现时通常配合位运算。

template <typename T>
T expow(T a, T b, const T P) {
    T ret = 1;
    for (; b; b >>= 1) {
        if (b & 1)
            ret = 1ll * ret * a % P;
        a = 1ll * a * a % P;
    }
    return ret;
}

从实质中可以猜想出矩阵也可以使用快速幂,然后发现这是对的。

注意到这个东西是可以直接套到矩阵上的,封装一下矩阵然后重载乘法运算符就可以了。

矩阵优化加速

你说的对,但是竞赛显然不会出板子矩阵题,那利用矩阵的性质我们可以干什么呢?

来看经典的 Fibonacci 数列问题,LG 上有类似的 P1939

我们都知道 Fibonacci 的递推式为:

$$ f_{i} = f_{i - 1} + f_{i - 2} $$

直接递推复杂度 $O(n)$,而且不滚动的话空间消耗巨大。我们考虑到 $f_{n}$ 只和 $f_{i - 1}, f_{i - 2}$ 有关,所以在递推的时候保留最近的两个 Fibonacci 数,就可以推导出下一个 Fibonacci 数了。

设 $F(n)$ 表示一个 $1 \times 2$ 的矩阵:$F(n) = \begin{bmatrix}f_{n} & f_{n + 1}\end{bmatrix}$。

把 $F(n - 1) = \begin{bmatrix}f_{n - 1} & f_{n}\end{bmatrix}$ 来计算出 $F(n)$,其实就是把上一个矩阵的第二个数搞到当前矩阵的第一个数来,把上一个矩阵第一、二列上的数字类加到当前矩阵的第二个数来,构造一个矩阵就可以了:

$$ F(n) = F(n - 1) \times \boldsymbol{A} \Leftrightarrow \begin{bmatrix}f_{n} & f_{n + 1}\end{bmatrix} = \begin{bmatrix}f_{n - 1} & f_{n}\end{bmatrix} \times \begin{bmatrix}0 & 1 \\ 1 & 1 \\ \end{bmatrix} $$

初始矩阵 $F(0) = \begin{bmatrix}0 & 1\end{bmatrix}$,目标 $F(n) = F(0) \times \boldsymbol{A}^{n}$,矩阵乘法满足结合律,快速幂来算,复杂度为 $O(2^3 \log n)$。

结论

归纳一下要点,对于一类满足:

  • 可抽象出 $n$ 长度的一维向量且向量在单位时间内发生一次变化
  • 变化的形式是一个线性的递推
  • 递推式一定,即在不同时间可能作用于不同数据但是式子本身保持不变
  • 向量变化时间,也即递推轮数很长,但向量长度 $n$ 不大

就可以考虑使用矩阵乘法来进行优化。这里有几个定义:

  • 长度为 $n$ 的一维向量:“状态矩阵”
  • 与状态矩阵相乘的一个固定的矩阵 $\boldsymbol{A}$:“转移矩阵”

状态矩阵中第 $x$ 个数对下一单位时间状态矩阵中的第 $y$ 个数产生影响,转移矩阵对应的 $(x, y)$ 位置就需要填充系数。

关键和难点在于“状态矩阵”,从递推式中搞出来之后,利用快速幂和矩阵乘法就可以写。

复杂度一般为 $O(n^3 \log T)$,其中 $T$ 为递推总轮数。

「馬鹿なわたしは歌うだけ」