矩阵加速

نُشر في عُدل في

矩阵加速是一种能将复杂度为$\Theta(n)$的递推式加速到$\Theta(k^3logn)$的方法(k为加速矩阵的大小)

1.原理

​ 假设我们有一满足递推关系的数列${a_i}$,已知其前$k$项(其中$k$就是加速矩阵的大小,具体为什么后面有解释)${a_1,a_2\cdots a_k}$

​ 要想求出$a_{k+1}$,我们可以构造出这样一个大小为$k\times k$的矩阵$A$使其满足
$$
\left[
\begin{matrix}
a_1 \
a_2 \
\vdots \
a_k
\end{matrix}
\right]
\times A

\left[
\begin{matrix}
a_2 \
a_3 \
\vdots \
a_{k+1}
\end{matrix}
\right]
$$
​ 再由矩阵乘法满足交换律可知想要求出$a_{k+n}$我们要把$A^n$乘上我们的原始矩阵就在矩阵的最后一行得到我们想要的$a_{k+n}$,即
$$
\left[
\begin{matrix}
a_1 \
a_2 \
\vdots \
a_k
\end{matrix}
\right]
\times A^n

\left[
\begin{matrix}
a_{n+1} \
a_{n+2} \
\vdots \
a_{n+k}
\end{matrix}
\right]
$$
很妙是吧(滑稽

​ 那么问题来了怎么求$A$呢

​ 由我上一篇博客

​ 矩阵快速幂

​ 可知如果我们把原始矩阵当作参数,所求矩阵当作方程的解,那么$A$就是我们的系数,那么我们不难得出
$$
\begin{cases}
aa_1+ba_2+\cdots +ca_k=a_2\
da_1+ea_2+\cdots +fa_k=a_3\
\vdots\
ga_1+ha_2+\cdots +ia_k=a_4
\end{cases}
\Rightarrow
A =
\left[\begin{matrix}
a&b\cdots c\
d&e\cdots f\
\vdots&\vdots \ddots \vdots \
g&h\cdots i
\end{matrix}\right]
$$
​ 现在我们有了矩阵$A$,可我们发现,这个$A$还要自乘$n$次,要知道,一次矩阵乘法可是$\Theta(k^3)$的复杂度,这时,矩阵快速幂再次发挥用处,我们通过他就可以把它降到$\Theta(logn)$综合起来我们就能得到文章开头的$\Theta(k^3logn)$的复 杂度

例子

​ 我们不如就拿最经典的递推式—–斐波那契数列来举例子吧

​ 我们知道斐氏数列的递推式是
$$
a_1=1,a_2=1\a_n=a_{n-1}+a_{n-2}
$$
​ 我们构造如下式子
$$
\left[
\begin{matrix}
1 , (a_1)\
1 , (a_2)\
\end{matrix}
\right]
\times A

\left[
\begin{matrix}
a_2 \
a_3
\end{matrix}
\right]
$$
​ 求$A$矩阵
$$
\begin{cases}
0 \times a_1 + 1 \times a_2 = a_2 \
1 \times a_1 + 1 \times a_2 = a_3
\end{cases}
$$
​ ∴$A$为
$$
\left[
\begin{matrix}
0&1\
1&1
\end{matrix}
\right]
$$
​ 所以我们如果要求$a_n$只需将$A$自乘$n-2$次再乘原来的矩阵就能得到$a_n$

​ 但是

​ 细心的读者可能会发现这个原始矩阵很特殊,我们也许不需要乘它就可以得到$a_n$我们试一下
$$
A^2=
\left[
\begin{matrix}
1&1\
1&2
\end{matrix}
\right]
\
A^3=
\left[
\begin{matrix}
1&2\
2&3
\end{matrix}
\right]
\
A^4=
\left[
\begin{matrix}
2&3\
3&5
\end{matrix}
\right]
$$
​ 发现什么了吗?

​ 我们不难看出当$A^{n}$时$(1,0),(0,1)$这两个数就是$a_{n}$!

​ 完整代码如下

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
#include<iostream>
using namespace std;
const long long mod = 1000000007;
struct mat
{
	long long m[5][5];
}a,e;

mat multi(mat a,mat b)
{
	mat c;
	for(int i=1;i<=2;i++)
		for(int j=1;j<=2;j++)
			c.m[i][j]=0;
	for(int i=1;i<=2;i++)
		for(int j=1;j<=2;j++)
			for(int k=1;k<=2;k++)
				c.m[i][j]=c.m[i][j]%mod+a.m[i][k]*b.m[k][j]%mod;
	return c;
}

mat ksm(mat x,long long p)
{
	mat res=e;
	while(p)
	{
		if(p&1)
			res=multi(res,a);
		a=multi(a,a);
		p>>=1;
	}
	return res;
}

long long n,m;

void init()
{
	a.m[1][1]=0; a.m[1][2]=1;
	a.m[2][1]=1; a.m[2][2]=1;
}

int main()
{
	for(int i=1;i<=2;i++)
		e.m[i][i]=1;
	scanf("%lld",&n);
	if(n==1||n==2)
	{
		cout<<1;
		return 0;
	}
	init();
	mat res=ksm(a,n);
	printf("%lld\n",res.m[1][2]%mod);
	return 0;
}

​