DAY_5（矩阵快速幂）

hycccccch

已于 2024-01-20 21:09:32 修改

阅读量404

点赞数 5

文章标签：矩阵算法

于 2024-01-20 01:32:03 首次发布

本文链接：https://blog.csdn.net/hycccccch/article/details/135709470

版权

本文详细介绍了矩阵快速幂算法，包括基本概念、实现原理、递归构造以及在多项式求和、区间和和二维递归方程中的应用。通过实例演示和代码展示，帮助读者理解并掌握这一高效计算技术。

摘要由CSDN通过智能技术生成

在讲矩阵快速幂之前，先来回顾快速幂：

快速幂

用于解决 a^n 当n很大时的情况。通常与取模同时应用。

a^n = a*a*a*a*.... 。难道我们真的要这么求吗？？？？？？

不可能的！！

经过快速幂算法，我们会成功地将时间复杂度从 O(n) 降低为 O(logn)

如何实现呢？举个栗子： 3*3*3*3*3*3*3*3 = (3*3)*(3*3)*(3*3)*(3*3) = ((3*3)*(3*3))*((3*3)*(3*3)) ,即 3^8 = 9^4 = 81^2 ,如果n为奇数， $a^n = a*a^{[n/2]}*a^{[n/2]}$ ,那么换成代码就是：

int fastpow(int m,int n,int k)
{
	int l=1;
	while(k!=0)
	{
		if(k%2)
		l*=m;
		m*=m;
		k/=2;
	}
	return l;
}

下面是矩阵快速幂！！！

矩阵快速幂

首先复习一下，矩阵乘法（洛谷是很好的！！）

（来源：洛谷P3390 【模板】矩阵快速幂题目背景）

矩阵乘法的代码：

for(ll i=1;i<=n;i++)
for(ll j=1;j<=n;j++)
for(ll k=1;k<=n;k++)
c[i][j]+=a[i][k]*b[k][j];

好了，有了矩阵乘法法则，我们就来看看对于有递推公式，但是数据范围很大的情况下如何解决问题：：

有如原始公式： f[n]=f[n-1]+f[n-2]

可以把它转化成矩阵相乘：： $\binom{f[n]}{f[n-1]}=\binom{1,1}{1,0}\binom{f[n-1]}{f[n-2]}$

我们这就完成了矩阵快速幂的构造（关于构造方法在下文）

所以我们现在就有了一个矩阵 $B=\binom{1,1}{1,0}$

一直递归下去，得到的最终结果是： $\binom{f[n]}{f[n-1]}=\binom{1,1}{1,0}^{n-1}\binom{f[1]}{f[0]}$

因此，想要得到 f[n] ,我们只需要求出 $B^{n-1}$ 就行

这就用到了矩阵快速幂

很简单，只需把快速幂变量换成定义矩阵的struct就行

juzhen fastpow(juzhen m,ll n,ll k)
{
	k--;
	juzhen l=m;
	while(k!=0)
	{
		if(k%2)
		l=cheng(l,m);//矩阵乘法
		m=cheng(m,m);//矩阵乘法
		k/=2;
	}
	return l;
}

接下来说说构造：

一、多项式

辅助矩阵一定是方阵，并且它的阶数等于递归方程的递归深度

$\bigl(\begin{smallmatrix} f[n]\\ f[n-1]\\ f[n-2] \end{smallmatrix}\bigr)=\bigl(\begin{smallmatrix} a & 0& b\\ 1 & 0 & 0\\ 0 & 1 & 0 \end{smallmatrix}\bigr)\bigl(\begin{smallmatrix} f[n-1]\\ f[n-2]\\ f[n-3] \end{smallmatrix}\bigr)$

二、有常数的多项式

由常数，当做多一阶处理

$\bigl(\begin{smallmatrix} f[n]\\ f[n-1]\\ f[n-2] \\c\end{smallmatrix}\bigr)=\bigl(\begin{smallmatrix} a & 0& b&1\\ 1 & 0 & 0&0\\ 0 & 1 & 0&0 \\0 &0&0&1\end{smallmatrix}\bigr)\bigl(\begin{smallmatrix} f[n-1]\\ f[n-2]\\ f[n-3]\\c \end{smallmatrix}\bigr)$

三、有变量的多项式

$f[n]=af[n-1]+bf(n-2)+n^{3}$

有变量，使用二次项定理将其分解，在当做多个阶来处理

因为 $n^{3}=(n-1+1)^{3}=(n-1)^{3}+3(n-1)^{2}+3(n-1)+1$

于是有了： $\bigl(\begin{smallmatrix} f[n]\\ f[n-1]\\ n^{3}\\ n^{2}\\ n\\1 \end{smallmatrix}\bigr)=\bigl(\begin{smallmatrix} 1&2 &1 &3 &3 &1 \\ 1& 0& 0& 0 & 0 & 0\\ 0 & 0 & 1 & 3 & 3 & 1\\ 0 & 0 & 0 & 1 & 2 & 1\\ 0 & 0 & 0 & 0 & 1 & 1\\ 0 & 0 & 0 & 0 & 0 & 1 \end{smallmatrix}\bigr)\bigl(\begin{smallmatrix} f[n-1]\\ f[n-2]\\ (n-1)^{3}\\ (n-1)^{2}\\ (n-1)\\ 1 \end{smallmatrix}\bigr)$

四、多项式的区间和

求 $f[a]+f[a+1]+f[a+2]+...+f[b-1]+f[b]$ 的值

（前缀和！！！）——想到啥了？用求和之差！！

设 $S[n]=\sum_{i=1}^{n}f[i]$ ，

则 $f[a]+f[a+1]+f[a+2]+...+f[b-1]+f[b]=S[b]-S[a-1]$

$S[n]=S[n-1]+f[n]$

故 $S[n]=S[n-1]+f[n-1]+f[n-2]+f[n-3]$

所以 $\bigl(\begin{smallmatrix} S[n]\\f[n]\\ f[n-1]\\ f[n-2] \end{smallmatrix}\bigr)=\bigl(\begin{smallmatrix} 1 & 1& 1&1\\ 1 & 0 & 0&0\\ 0 & 1 & 0&0 \\0 &0&1&0\end{smallmatrix}\bigr)\bigl(\begin{smallmatrix} S[n-1]\\f[n-1]\\ f[n-2]\\ f[n-3] \end{smallmatrix}\bigr)$

五、含乘法项

f[n]=x*f[n-1]+y*f[n-2] ，求 $S[n]=\sum f[n]^{2}$

推导：

$S[n]=S[n-1]+f[n]^{2}$ ， $f[n]^{2}=x^{2}f[n-1]^{2}+2xyf[n-1]f[n-2]+y^2f[n-2]^2$

f[n]f[n-1]=x*f[n-1]^2+y*f[n-1]f[n-2]

所以 $\bigl(\begin{smallmatrix} S[n]\\f[n]^2\\ f[n]f[n-1]\\ f[n-1]^2 \end{smallmatrix}\bigr)=\bigl(\begin{smallmatrix} 1 & x^2& 2xy&y^2\\ 0 & x^2 & 2xy&y^2\\ 0 & x & y&0 \\0 &1&0&0\end{smallmatrix}\bigr)\bigl(\begin{smallmatrix} S[n-1]\\f[n-1]^2\\f[n-1]f[n-2]\\ f[n-2]^2\end{smallmatrix}\bigr)$

六、二维递归方程

$f[n,l]=f[n-1,l-1]+f[n-1,l]$

对于二维递归方程，

构造：： $\bigl(\begin{smallmatrix} f[n,1]\\ f[n,2]\\ ...\\ f[n,l] \end{smallmatrix}\bigr)=\bigl(\begin{smallmatrix} 1& 0 &... & 0 & 1\\ 1& 1 & ... & 0 & 0\\ ...& ... & ... & 1 & 0\\ 0& 0 & ... & 1 & 1 \end{smallmatrix}\bigr)\bigl(\begin{smallmatrix} f[n-1,1]\\ f[n-1,2]\\ ...\\ f[n-1,l] \end{smallmatrix}\bigr)$

下面是一些题目：：

P3390 【模板】矩阵快速幂

板子题，主要就是快速幂和矩阵乘法两个函数，没有太大思维量

代码：

#include<bits/stdc++.h>
using namespace std;
#define ll long long
const int N=105;
const int md=1e9+7; 
ll n,k;
struct juzhen
{
	ll a[N][N];
}p;
juzhen cheng(juzhen s,juzhen b)
{
	juzhen c;
	memset(c.a,0,sizeof(c.a));
	for(ll i=1;i<=n;i++)
	for(ll j=1;j<=n;j++)
	for(ll k=1;k<=n;k++)
	c.a[i][j]=(c.a[i][j]%md+(s.a[i][k]*b.a[k][j])%md)%md;
	return c;
}
juzhen fastpow(juzhen m,ll n,ll k)
{
	k--;
	juzhen l=m;
	while(k!=0)
	{
		if(k%2)
		l=cheng(l,m);
		m=cheng(m,m);
		k/=2;
	}
	return l;
}
int main()
{
	cin>>n>>k;
	for(ll i=1;i<=n;i++)
	for(ll j=1;j<=n;j++)
	{
		cin>>p.a[i][j];
		if(k==0)
		{
			p.a[i][j]=0;
			p.a[i][i]=1;
		}
	}
	juzhen ans=fastpow(p,n,k);
	for(ll i=1;i<=n;i++)
	{
		for(ll j=1;j<=n;j++)
		cout<<ans.a[i][j]%md<<" ";
		cout<<endl;
	}
}

P1939 矩阵加速（数列）

根据题目给定的递推公式，先构造一个矩阵，然后进行矩阵幂运算，最后求an也就是最后矩阵第一行的前两个数相加即可（一开始以为是第一个和第三个，wa了好几次。。）

代码：

#include<bits/stdc++.h>
using namespace std;
#define ll long long
const int md=1e9+7;
ll t;
ll n;
struct jz
{
	ll a[4][4];
}h;
jz quick_pow(jz p,jz c)
{
	jz b;
	memset(b.a,0,sizeof(b.a));
	for(ll i=1;i<=3;i++)
	for(ll j=1;j<=3;j++)
	for(ll k=1;k<=3;k++)
	b.a[i][j]=(b.a[i][j]%md+(p.a[i][k]*c.a[k][j])%md)%md;
	return b;
}
ll cf(jz q,ll n)
{
	if(n==1||n==2||n==3)
	return 1;
	ll k=n-3;
	jz o=h;
	while(k)
	{
		if(k%2)
		q=quick_pow(q,o);
		o=quick_pow(o,o);
		k/=2;
	}
	return q.a[1][1]+q.a[1][2];
}
int main()
{
	cin>>t;
	memset(h.a,0,sizeof(h.a));
	h.a[1][1]=1;
	h.a[1][3]=1;
	h.a[2][1]=1;
	h.a[3][2]=1;
	for(int i=1;i<=t;i++)
	{
		cin>>n;
		cout<<cf(h,n)%md<<endl;
	}
}

P1397 [NOI2013] 矩阵游戏

其实这题不是特别难，放在蓝题里面属实有点水分，不过数据范围，额，也的确恶心的可以

题目给了三个递推公式，可是根据这个构造不出矩阵来，于是需要进行推导：

所以一开始可以得到 $\bigl(\begin{smallmatrix} f[n,m]\\ b \end{smallmatrix}\bigr)=\bigl(\begin{smallmatrix} a & 1\\ 1 &0 \end{smallmatrix}\bigr)\bigl(\begin{smallmatrix} f[n,m-1]\\ b \end{smallmatrix}\bigr)$

然后一直往下递推，就有： $\bigl(\begin{smallmatrix} f[n,m]\\ b \end{smallmatrix}\bigr)=\bigl(\begin{smallmatrix} a & 1\\ 1 &0 \end{smallmatrix}\bigr)^{m-1}\bigl(\begin{smallmatrix} f[n,1]\\ b \end{smallmatrix}\bigr)$

可以设 $\bigl(\begin{smallmatrix} a & 1\\ 1 &0 \end{smallmatrix}\bigr)^{m-1}=B=\bigl(\begin{smallmatrix} x1 & y1\\ ... &0 \end{smallmatrix}\bigr)$

那么 $f[n,m]=x1*f[n,1]+y1*b$

再根据第三个递推式得： $f[n,1]=c*f[n-1,m]+d$

代入上面第二行的式子： $f[n,m]=c*x1*f[n-1,m]+(x1*d+y1*b)$

所以现在就可以来构造矩阵了： $\bigl(\begin{smallmatrix} f[n,m]\\ x1*d+y1*b \end{smallmatrix}\bigr)=\bigl(\begin{smallmatrix} c*x1 & 1\\ 1 &0 \end{smallmatrix}\bigr)^{n-1}\bigl(\begin{smallmatrix} f[1,m]\\ x1*d+y1*b \end{smallmatrix}\bigr)$

再来设 $\bigl(\begin{smallmatrix} c*x1 & 1\\ 1 &0 \end{smallmatrix}\bigr)^{n-1}=C=\bigl(\begin{smallmatrix} x2 & y2\\ ...& 0 \end{smallmatrix}\bigr)$

所以 $f[n,m]=x2*f[1,m]+y2*(x1*d+y1*b)$

又因为 $f[1,m]=x1*f[1,1]+y1*b=x1+y1*b(f[1,1]=1)$

最后代入得到 $f[n,m]=x1*x2+x2*y1*b+x1*y2*d+y1*y2*b$

代码：

#include<bits/stdc++.h>
using namespace std;
#define ll long long 
const int md=1e9+7;
const int mdd=md-1;
ll n,m,a,b,c,d;
struct jz
{
	ll matrix1[3][3];
	ll matrix2[3][3];
	ll matrix3[3][3];
}p;
inline void readBigInt(ll &ans1, ll &ans2)
{
	ans1=0,ans2=0;
	char ch=getchar();
	while(!isdigit(ch))
	ch=getchar();
	while(isdigit(ch))
	{
		ans1=(ans1<<3)+(ans1<<1)+(ch^48); // 返回两种取模结果
		ans1%=md;
		ans2=(ans2<<3)+(ans2<<1)+(ch^48);
		ans2%=mdd;
		ch=getchar();
	}
}
jz quick_pow1(jz p,jz c)
{
	jz x;
	memset(x.matrix1,0,sizeof(x.matrix1));
	for(ll i=1;i<=2;i++)
	for(ll j=1;j<=2;j++)
	for(ll k=1;k<=2;k++)
	x.matrix1[i][j]=(x.matrix1[i][j]%md+(p.matrix1[i][k]*c.matrix1[k][j])%md)%md;
	return x;
}
jz quick_pow3(jz p,jz c)
{
	jz x;
	memset(x.matrix3,0,sizeof(x.matrix3));
	for(ll i=1;i<=2;i++)
	for(ll j=1;j<=2;j++)
	for(ll k=1;k<=2;k++)
	x.matrix3[i][j]=(x.matrix3[i][j]%md+(p.matrix3[i][k]*c.matrix3[k][j])%md)%md;
	return x;
}
jz cf1(jz q,ll k)
{
	jz h=q;
	while(k)
	{
		if(k%2)
		h=quick_pow1(h,q);
		q=quick_pow1(q,q);
		k/=2;
	}
	return h;
}
jz cf3(jz q,ll k)
{
	jz h=q;
	while(k)
	{
		if(k%2)
		h=quick_pow3(h,q);
		q=quick_pow3(q,q);
		k/=2;
	}
	return h;
}
int main()
{
	ll n1,n2,m1,m2;
	readBigInt(n1,n2),readBigInt(m1,m2);
	cin>>a>>b>>c>>d;
	a==1? (n=n1,m=m1):(n=n2,m=m2); // 特判
	p.matrix1[1][1]=a%md;
	p.matrix1[1][2]=1;
	p.matrix1[2][2]=1;
	jz ans1=cf1(p,m-2);
	ll x1=ans1.matrix1[1][1]%md;
	ll y1=ans1.matrix1[1][2]%md;
	p.matrix3[1][1]=(x1*c)%md;
	p.matrix3[1][2]=1;
	p.matrix3[2][2]=1;
	jz ans3=cf3(p,n-2);
	ll x2=ans3.matrix3[1][1]%md;
	ll y2=ans3.matrix3[1][2]%md;
	ll cs=((x1*d)%md+(y1*b)%md)%md;
	cout<<((((x1+(y1*b)%md)%md)*x2)%md+(y2*cs)%md)%md<<endl;
}

哦对，一开始只过了一半数据点，原因是忽略了数据范围。。。

。。。。。。。。。。。。。。这太抽象了！！！！！！！

然后题解里神犇们用了字符串读入 并取模！！！总之 t q l ，然后CV到了自己代码里嘿嘿~~

最后hack了两个数据点：

34578734657863487 38465873465876348 1 23734637 3892734 3849

467549493

1145141919810 1000000007 2 1 1 2

283823589

不知道为啥，请求大佬指正！！

~~完结撒花~~~

hycccccch

关注

5
点赞
踩
7

收藏

觉得还不错? 一键收藏
0
评论
DAY_5（矩阵快速幂）

（来源：洛谷P3390 【模板】矩阵快速幂题目背景）很简单，只需把快速幂变量换成定义矩阵的struct就行。我们这就完成了矩阵快速幂的构造（关于构造方法在下文）经过快速幂算法，我们会成功地将时间复杂度从。首先复习一下，矩阵乘法（洛谷是很好的！难道我们真的要这么求吗？将其分解，在当做多个阶来处理。所以我们现在就有了一个矩阵。未完待续......
复制链接

扫一扫