矩阵快速幂

矩阵快速幂

前言

第一次听说矩阵快速幂是在NOI 2023 冬令营的课堂上，当时自己以为很高深，有关线代就一只没学，没了解。直到有一天……

引入

动机

下图，这是我在学习概率DP时，看aqx的PPT中的一道例题，这里咱们暂且先不讨论本题的做法，我只要是想通过这个题来引入主题，看到本题的Task 2，是不是很难想啊……我也是问过wvr之后才知道正解是矩阵快速幂，所以下定决心今天共同学习一下矩阵快速幂

当然，我们不能通过上面的例子来引入今天的主题，虽然上面的例子是我学习矩阵快速幂的动机

那么咱们来看下面的这个例子

P3390 【模板】矩阵快速幂

很简单的一道题哈

那么额……通过这道题我们得出了一个什么事情？

见下

定义

矩阵快速幂就是将快速幂运用到矩阵乘法上的算法

其最常见的应用是加速矩阵的计算，以较低的时间复杂度内得到递推公式的第 \(n\) 项

前置知识

矩阵乘法

矩阵 \(A_{n \times m}\) 本质上是一个 \(n\) 行 \(m\) 列的二维数组

对于矩阵乘法 \(C_{n \times s} = A_{n \times m} \times B_{m \times s}\) 来说，\(C_{i,j} = \Sigma A_{i,k} \times B_{k,j} \space \space (1 \le k \le m)\)

通俗地来讲，对于矩阵 \(C\) 第一行第一列的元素，是矩阵 \(A\) 的第一行元素，与矩阵 \(B\) 的第一列元素，对应相乘再相加的和的结果……然后类推下去

代码如下：

struct matrix {
    ll x[105][105];
    matrix () {memset(x, 0, sizeof(x));}
};
matrix multiply(matrix &a, matrix &b) {
    matrix c;
    for (int i = 1; i <= n; i++)
        for (int j = 1; j <= n; j++)
            for (int k = 1; k <= n; k++)
                c.x[i][j] += a.x[i][k] * b.x[k][j];
    return c;
}

快速幂

快速幂（Exponentiation by squaring，平方求幂）是一种简单而有效的小算法

它可以 \(\Theta (\log n)\) 的时间复杂度计算乘方

让我们来回忆一下快速幂的模板

while (n > 0) {
    if (n & 1) res *= a;
    a = a * a;
    n >>= 1;
}

单位矩阵

定义：

一个主对角线（左上-右下）元素均是1，其余元素均是0

记作符号 \(E\) 或 \(I\)

相对应的 \(n\) 维单位矩阵可记作 \(E_n\) 或 \(I_n\)

例如下图就是一个标准的单位矩阵

\[\begin{bmatrix} 1&0&0&0\\ 0&1&0&0\\ 0&0&1&0\\ 0&0&0&1 \end{bmatrix} \]

应用详解

通过上面的介绍，有一件事情变得非常显然，那就是矩阵快速幂 = 矩阵乘法 + 快速幂

矩阵快速幂只需要将快速幂里的数字放进矩阵即可

注：我们需要首先构造一个单位矩阵（任何矩阵乘单位矩阵结果都是本身，相当于数字1）在进行幂的运算

代码如下：

matrix mpow (matrix &a, ll m) { 				//矩阵a的m次方 
    matrix res;
    for (int i = 1; i <= n; i++) res.x[i][i] = 1; //单位矩阵
    while (m > 0) {
        if (m & 1) res = multiply(res, a);
        a = multiply(a, a);
        m >>= 1;
    }
    return res;
}

另外，大部分题目会在求矩阵的幂的同时，要求模一个数，一般我们选择在矩阵乘法的同时，进行取模操作

原则

记住，原理是对于递推式有优化哈，比如你写DP的转移方程，是吧？

斐波那契数列

对于斐波那契数列，我们在小学二年级就学过的

\(f[0] = 0, \space f[1] = 1, \space f[2] = 2 \space ……\space f[n] = f[n - 1] + f[n - 2]\)

正常的线性递推，时间复杂度是 \(\Theta (n)\) 的，当 \(n\) 较大时，难以计算

而使用矩阵快速幂之后，时间复杂度可以达到 \(\Theta (\log n)\)

那么具体如何操作？见下：

我们需要构造俩个矩阵 \(A, \space B\)

对于矩阵 \(A\) 必须含有俩个初始值（这里可以考虑因为 \(f[n] = f[n - 1] + f[n - 2]\)）

即是：\(f[0]\) 和 \(f[1]\) 的值

所以

\[A = \begin{bmatrix} 0\\ 1\\ \end {bmatrix} \]

对于矩阵 \(B\) 我们并不知道是个啥子？

所以可以先设出来：

\[B = \begin{bmatrix} a & b\\ c & d\\ \end{bmatrix} \]

然后给 \(A\) 左乘以 \(B\) 得到了，\(f[0] + f[1] = f[2]\) ，就是从 \(0\) 和 \(1\) 推出 \(2\) 是多少

矩阵表达为：

\[\begin{bmatrix} a & b\\ c & d\\ \end{bmatrix} \times \begin{bmatrix} 0\\ 1\\ \end{bmatrix} = \begin{bmatrix} 1\\ 1\\ \end{bmatrix} \]

重复第一次的计算，从 \(f[1] + f[2] = f[3]\) ，就是从 \(1\) 和 \(2\) 推出 \(3\) 是多少的过程

即给新矩阵左乘以 \(B\) 得到了：

\[\begin{bmatrix} a & b\\ c & d\\ \end{bmatrix} \times \begin{bmatrix} 1\\ 1\\ \end{bmatrix} = \begin{bmatrix} 1\\ 2\\ \end{bmatrix} \]

现在我们得到了俩个矩阵表达式子，这样可以解得 \(a,b,c,d\) 分别是多少

最终得到：

\[B = \begin{bmatrix} 0 & 1\\ 1 & 1\\ \end{bmatrix} \]

综上，我们已经推出来了矩阵 \(B\) ，接下来就是一次一次的往后乘，求第 \(m\) 项就是求 \(A \times B^m\) ，乘的过程就用矩阵乘法

总结

使用矩阵快速幂的时候，一般构造俩个矩阵 \(A, \space B\) ，满足一下条件：

1. 第一个矩阵 \(A\) 含且只含有初始值
2. 第一个矩阵左乘第二个矩阵，得到的新矩阵的大小不变
3. 每一次左乘或者说是每一次特定的线性计算后，第 \(i\) 项都会变成第 \(i + 1\) 项

写在最后

我终于写完了啊，写博客有点耗时间啊，这篇估计写了 \(3-4\) 个小时，主要是 \(Markdown\) 不太熟悉，还有就是过程中边学边写，有点慢了，不过问题不大，最终的效果还是很好的。。。