从fibonacci数列开始说起
给定一个nn,求fibonacci数列的第nn项模一个数PP,即F_n \mod PFnmodP。
题目:洛谷1962
较朴素的方法:递推求数
int fibonacci(int n) {
static int f[N];
f[1] = 1;
for (int i = 2; i <= n; ++i)
f[i] = f[i - 1] + f[i - 2] % P;
return f[n];
}
这样,我们可以发现其时间复杂度与空间复杂度都是O(n)O(n),当然,我们可以发现每个数都只需要它之前的两个数,所以我们可以通过这一性质进行空间压缩,将空间复杂度降到O(1)O(1)。这个效率确实不错,但是当n的取值范围提升到10^{12}1012这个数量级时便无能为力了。
事实上,有一个方法可以以O(\log n)O(logn)的效率将F_nFn算出。
将数列递推式转化为线性变换
我们可以发现,将空间压缩后的每个时刻状态,可以当做一个矢量:
A\_n = \left( \begin{array}{cc} F\_{n-1} & F\_n \end{array} \right)A_n=(F_n−1F_n)
在转化成下一个状态的过程中:
A\_{n+1} = \left( \begin{array}{cc} F\_{n} & F\_{n+1} \end{array} \right) = \left( \begin{array}{cc} A\_{n2} & A\_{n1}+A\_{n2} \end{array} \right)A_n+1=(F_nF_n+1)=(A_n2A_n1+A_n2)
我们便可以将其改写成一个矩阵变换,设
M = \left( \begin{array}{cc} 0 & 1 \\ 1 & 1 \end{array} \right)M=(0111)
则有
A_{n+1} = A_n MAn+1=AnM
根据归纳法,我们可以得知
A_n=A_1 M^{n-1}An=A1Mn−1
进而有
F\_n=(A\_1 M^{n-1})\_2=(M^{n-1})\_{2,2}F_n=(A_1Mn−1)_2=(Mn−1)_2,2
那我们要问了,矩阵乘法的复杂度是O(r^3)O(r3)的,这里rr为2即乘方运算大致是8n8n次加法,而原本的求法是nn次加法,这岂不是常数反而更大了?实际上,矩阵乘方和数字的普通乘方运算,都可以通过快速幂来进行加速。那么便可以通过O(\log n)O(logn)的复杂度完成乘方运算。接下来,我将给出证明。
代码样例:
#include <cstdio>
#include <cstring>
using namespace std;
typedef long long ll;
struct matrix {
int a[2][2];
matrix()
{ memset(a, 0, sizeof(a)); }
matrix(const matrix& x)
{ memcpy(a, x.a, sizeof(a)); }
matrix operator*(const matrix& x) const;
};
const int P = 1e9+7;
int main() {
matrix ans, m;
ll n;
scanf("%lld", &n);
--n;
ans.a[0][0] = ans.a[1][1] = 1;
m.a[0][1] = m.a[1][0] = m.a[1][1] = 1;
while (n) {
if (n & 1)
ans = ans * m;
m = m * m;
n >>= 1;
}
printf("%d\n", ans.a[1][1]);
return 0;
}
matrix matrix::operator*(const matrix& x) const {
matrix ret;
for (int i = 0; i < 2; ++i)
for (int j = 0; j < 2; ++j)
for (int k = 0; k < 2; ++k)
ret.a[i][j] = (((ll) a[i][k] * x.a[k][j]) + ret.a[i][j]) % P;
return ret;
}
矩阵乘法的基本定理
矩阵乘法的结合律
设有矩阵A,B,CA,B,C分别的大小为n \times m, m \times p, p \times qn×m,m×p,p×q。求证(AB)C=A(BC)(AB)C=A(BC),进一步为(AB)C_{i,j}=A(BC)_{i,j}(AB)Ci,j=A(BC)i,j
根据矩阵乘法的定义,有
AB\_{i,j} = \sum\_{k=1}^m A\_{i,k} \times B\_{k,j}AB_i,j=∑_k=1mA_i,k×B_k,j
故
$$ ```cpp \begin{eqnarray} (AB)C_{i,j} & = & \sum_{l=1}^p AB_{i,l} \times C_{l,j} \\ & = & \sum_{l=1}^p (\sum_{k=1}^m A_{i,k} \times B_{k,l}) \times C_{l,j} \\ & = & \sum_{l=1}^p \sum_{k=1}^m A_{i,k} \times B_{k,l} \times C_{l,j} \\ & = & \sum_{k=1}^m \sum_{l=1}^p A_{i,k} \times B_{k,l} \times C_{l,j} \\ & = & \sum_{k=1}^m (\sum_{l=1}^p B_{k,l} \times C_{l,j}) \times A_{i,k} \\ & = & \sum_{k=1}^m A_{i,k} \times BC_{k,j} \\ & = & A(BC)_{i,j} \end{eqnarray} ``` $$
证毕。
矩阵快速幂的正确性
设方阵AA(行数为rr)和自然数nn,可以将nn分解为一个mm位的二进制数
\overline{b_mb_{m-1}...b_1b_0}bmbm−1...b1b0,则
$$ \begin{eqnarray} A^n ```cpp & = & A^{\sum_{i=0}^m b_i 2^i} \\ & = & \prod_{i=0}^m A^{b_i 2^i} \\ & = & \prod_{b_i = 1} A^{2^i} \end{eqnarray} ``` $$
从第一步到第二步的变换是不平凡的,它只有基于结合律才成立。因此,我们可以发现通过反复平方法可以以O(r^3 \log n)O(r3logn)的效率求出A^nAn。
重定义运算符
矩阵乘法中的运算涉及两种操作,即乘法和加法,我们在对结合律进行证明的过程中可以发现,我们依赖了:两种运算符各自符合交换律、结合律,乘法对加法符合分配率。形式化描述,即有运算\oplus⊕加法与\otimes⊗乘法,只要它们符合以下性质
$$ ```cpp \begin{eqnarray} x \oplus y & = & y \oplus x \\ (x \oplus y) \oplus z & = & x \oplus (y \oplus z) \\ x \otimes y & = & y \otimes x \\ (x \otimes y) \otimes z & = & x \otimes (y \otimes z) \\ x \otimes (y \oplus z) & = & (x \otimes y) \oplus (x \otimes z) \end{eqnarray} ``` $$
那么将该运算替换掉原本矩阵乘法的加法和乘法,交换律依然成立,快速幂也仍然可以使用。比如说,我们知道带余加法和带余乘法仍然满足该性质,所以我们也可以置\times_p \rightarrow \otimes, +_p \rightarrow \oplus×p→⊗,+p→⊕。
当然,如果还想要让这种重定义的矩阵可逆等,就需要该矩阵的元素所属集合与两种运算构成域(field),也就是还要有逆元和单位元。
矩阵快速幂解决最短路问题
我们尝试置+ \rightarrow \otimes, \min \rightarrow \oplus+→⊗,min→⊕后,可以发现这仍然符合矩阵运算所需要的性质,而这时
AB\_{i,j} = \bigoplus\_{k=1}^n A\_{i,k} \otimes B\_{k,j} = \min\_{k=1}^n (A\_{i,k} + B\_{k,j})AB_i,j=⨁_k=1nA_i,k⊗B_k,j=min_k=1n(A_i,k+B_k,j)
(这个minmin的特别用法是我自己编的,嘻嘻)
我们惊奇地发现,这就等同于一次边(i,j)(i,j)的松弛操作!于是我们对图的邻接矩阵AdjAdj计算n-1n−1次幂,就可以得出两两边的最短路了。其中如果(u,v)(u,v)无边,则置Adj_{u,v} \leftarrow +\inftyAdju,v←+∞。计算Adj^{n-1}Adjn−1的时间复杂度为O(n^3 \log n)O(n3logn),虽然逊色于floyd算法的O(n^3)O(n3),但是这也是一个非常启发式的方法。
【推荐】国内首个AI IDE,深度理解中文开发场景,立即下载体验Trae
【推荐】编程新体验,更懂你的AI,立即体验豆包MarsCode编程助手
【推荐】抖音旗下AI助手豆包,你的智能百科全书,全免费不限次数
【推荐】轻量又高性能的 SSH 工具 IShell:AI 加持,快人一步
· 单线程的Redis速度为什么快?
· 展开说说关于C#中ORM框架的用法!
· Pantheons:用 TypeScript 打造主流大模型对话的一站式集成库
· SQL Server 2025 AI相关能力初探
· 为什么 退出登录 或 修改密码 无法使 token 失效