洛谷-矩阵快速幂-(矩阵加速整理)

可爱美少女

已于 2022-09-15 20:49:16 修改

阅读量326

点赞数 1

分类专栏：数论动态规划文章标签：矩阵算法

于 2022-09-13 17:21:16 首次发布

本文链接：https://blog.csdn.net/m0_52398496/article/details/126836798

版权

动态规划同时被 2 个专栏收录

84 篇文章 2 订阅

订阅专栏

数论

60 篇文章 1 订阅

订阅专栏

矩阵快速幂

题意：
就是给你一个n×n的矩阵，问你n的k次幂是多少。

思考：

就是只要明白快速幂的本质是什么，这个就很好办。比如求a^b，这个b完全可以由二进制表达出来，比如11 = 1011，所以求a^b就可以分成多个a^{2的多少次幂}的乘积。比如3¹¹= 3⁸⁺²⁺¹ = 3⁸×3²×3¹。所以对b进行求每一位是否为1的同时，a一直都平方，如果b的这一位为1，那么正好sum*当前的a。
所以矩阵和这个一样的，对于求矩阵快速幂的时候，肯定是行列相同的，要不然没法×。对于两个矩阵分别为(a行,b列)，(c行,d列)，那么在保证b=c情况下，会得到一个(a行,d列)的矩阵。
在结构体里写东西的时候，不能直接初始化任何东西，都必须要在函数里去初始化。当然在结构体内部写函数也可以调用结构体里面的函数。反正就是只要不是定义变量，那么别的都要在函数里写。

代码：

#include<bits/stdc++.h>
#define fi first
#define se second
#define pb push_back
#define db double
#define int long long
#define PII pair<int,int >
#define mem(a,b) memset(a,b,sizeof(a))
#define IOS std::ios::sync_with_stdio(false),cin.tie(0),cout.tie(0);

using namespace std;
const int mod = 1e9+7,inf = 1e18;
const int N = 2e5+10,M = 2010;

int T,n,m,k;
int va[N];

struct matx{
	int n,m,a[M][M],b[M][M],c[M][M];
	void init()
	{
		for(int i=1;i<=n;i++)
		for(int j=1;j<=n;j++) c[i][j] = (i==j);
	}
	void ksm1()
	{
		for(int i=1;i<=n;i++)
		{
			for(int j=1;j<=n;j++)
			b[i][j] = c[i][j];
		}
		for(int i=1;i<=n;i++)
		{
			for(int j=1;j<=n;j++)
			{
				c[i][j] = 0;
				for(int k=1;k<=n;k++)
				c[i][j] = (c[i][j]+a[i][k]*b[k][j]%mod)%mod;
			}
		}
	}
	void ksm2()
	{
		for(int i=1;i<=n;i++)
		{
			for(int j=1;j<=n;j++)
			b[i][j] = a[i][j];
		}
		for(int i=1;i<=n;i++)
		{
			for(int j=1;j<=n;j++)
			{
				a[i][j] = 0;
				for(int k=1;k<=n;k++)
				a[i][j] = (a[i][j]+b[i][k]*b[k][j]%mod)%mod;
			}
		}
	}
	void ksm()
	{
		while(m)
		{
			if(m&1) ksm1();
			ksm2();m >>= 1;
		}
	}
}mat;


signed main()
{
	IOS;
	cin>>n>>m;
	for(int i=1;i<=n;i++)
	{
		for(int j=1;j<=n;j++)
		cin>>mat.a[i][j];
	}
	mat.n = n,mat.m = m;mat.init();mat.ksm();
	for(int i=1;i<=n;i++)
	{
		for(int j=1;j<=n;j++)
		cout<<mat.c[i][j]<<" ";
		cout<<"\n";
	}
	return 0;
}

矩阵加速（数列）

题意：
给你一个数列a，1<=i<=3时a[i] = 1,i>=4时a[i] = a[i-1]+a[i-3]。现在让你求这个数列的第n项，问你是多少。

思考：

既然给你了3三项，那么就要构造一个目标矩阵，那么可以构造为(F[i],F[i-1],F[i-2])。那么这个矩阵如何从前面转移呢？根据题目给的式子可以推出：
f[i]=f[i−1]×1+f[i−2]×0+f[i−3]×1
f[i−1]=f[i−1]×1+f[i−2]×0+f[i−3]×0
f[i−2]=f[i−1]×0+f[i−2]×1+f[i−3]×0
所以可以构造出转移矩阵：
(1,0,1)
(1,0,0)
(0,1,0)
那么(F[i],F[i-1],F[i-2])*转移矩阵 = (F[i+1],F[i],F[i-1])。
所以从第4项开始，就要从第3项，左×转移矩阵，求到第n项，那么就要左乘转移矩阵的n-3次方。所以求个矩阵快速幂后再乘以(F[3],F[2],F[1])就可以了。

代码：

#include<bits/stdc++.h>
#define fi first
#define se second
#define pb push_back
#define db double
#define int long long
#define PII pair<int,int >
#define mem(a,b) memset(a,b,sizeof(a))
#define IOS std::ios::sync_with_stdio(false),cin.tie(0),cout.tie(0);

using namespace std;
const int mod = 1e9+7,inf = 1e18;
const int N = 2e5+10,M = 10;

int T,n,m,k;
int va[M][M];
int anw[M][M];

struct matx{
	int n,m,a[M][M],b[M][M],c[M][M];
	void init()
	{
		for(int i=1;i<=n;i++)
		for(int j=1;j<=n;j++) c[i][j] = (i==j);
	}
	void ksm1()
	{
		for(int i=1;i<=n;i++)
		{
			for(int j=1;j<=n;j++)
			b[i][j] = c[i][j];
		}
		for(int i=1;i<=n;i++)
		{
			for(int j=1;j<=n;j++)
			{
				c[i][j] = 0;
				for(int k=1;k<=n;k++)
				c[i][j] = (c[i][j]+a[i][k]*b[k][j]%mod)%mod;
			}
		}
	}
	void ksm2()
	{
		for(int i=1;i<=n;i++)
		{
			for(int j=1;j<=n;j++)
			b[i][j] = a[i][j];
		}
		for(int i=1;i<=n;i++)
		{
			for(int j=1;j<=n;j++)
			{
				a[i][j] = 0;
				for(int k=1;k<=n;k++)
				a[i][j] = (a[i][j]+b[i][k]*b[k][j]%mod)%mod;
			}
		}
	}
	void ksm()
	{
		while(m)
		{
			if(m&1) ksm1();
			ksm2();m >>= 1;
		}
	}
}mat;

signed main()
{
	IOS;
	cin>>T;
	while(T--)
	{
		cin>>n;
		if(n<=3) cout<<1<<"\n";
		else
		{
			mat.a[1][1] = 1,mat.a[1][2] = 0,mat.a[1][3] = 1;
			mat.a[2][1] = 1,mat.a[2][2] = 0,mat.a[2][3] = 0;
			mat.a[3][1] = 0,mat.a[3][2] = 1,mat.a[3][3] = 0;
			mat.n = 3,mat.m = n-3;mat.init();mat.ksm();
			va[1][1] = 1,va[2][1] = 1,va[3][1] = 1;
			for(int i=1;i<=3;i++)
			{
				for(int j=1;j<=1;j++)
				{
					anw[i][j] = 0;
					for(int k=1;k<=3;k++)
					anw[i][j] = (anw[i][j]+mat.c[i][k]*va[k][j]%mod)%mod;
				}
			}
			cout<<anw[1][1]<<"\n";
		}
	}
	return 0;
}

题意：
就是给你n个点，刚开始在第1点，然后问你经过k天后也回到1点，中间可以随便走，但是不能有两个相邻天并且点也相同的情况。问你方案数有多少，并取模。

思考：

刚开始想了想推了推式子推不出来，然后就想了想dp，这个呢dp就是要两种状态dp[i][0]代表到第i天并且不在1号点，dp[i][1]代表在第i天在第一号点，那么转移就是dp[i][0] = dp[i-1][1]×(n-1)+dp[i-1][0]×(n-2)，dp[i][1] = dp[i-1][0];
但是k很大，这跑一遍dp直接超时了。不过看递推式子，无非就是从前面的来转移，那么就可以用矩阵加速。
初始矩阵：(F[1][0],F[1][1]) = (0,1)
目标矩阵：(F[k][0],F[k][1])。
转移矩阵：
(n-2,n-1)
(1,0)
由于题目说k天后，所以矩阵是k次幂，并不是问你第k天是多少，如果是第k天，那么次幂就是k-1次方。一般初始矩阵都是刚开始可以看出来的手推的，目标矩阵就是dp的定义，主要就是如何根据转移方程推出转移矩阵。最后求出来次幂乘以初始矩阵就可以了。

代码：

#include<bits/stdc++.h>
#define fi first
#define se second
#define pb push_back
#define db double
#define int long long
#define PII pair<int,int >
#define mem(a,b) memset(a,b,sizeof(a))
#define IOS std::ios::sync_with_stdio(false),cin.tie(0),cout.tie(0);

using namespace std;
const int mod = 998244353,inf = 1e18;
const int N = 2e5+10,M = 10;

int T,n,m,k;
int va[M][M];
int anw[M][M];

struct matx{
	int n,m,a[M][M],b[M][M],c[M][M];
	void init()
	{
		for(int i=1;i<=n;i++)
		for(int j=1;j<=n;j++) c[i][j] = (i==j);
	}
	void ksm1()
	{
		for(int i=1;i<=n;i++)
		{
			for(int j=1;j<=n;j++)
			b[i][j] = c[i][j];
		}
		for(int i=1;i<=n;i++)
		{
			for(int j=1;j<=n;j++)
			{
				c[i][j] = 0;
				for(int k=1;k<=n;k++)
				c[i][j] = (c[i][j]+a[i][k]*b[k][j]%mod)%mod;
			}
		}
	}
	void ksm2()
	{
		for(int i=1;i<=n;i++)
		{
			for(int j=1;j<=n;j++)
			b[i][j] = a[i][j];
		}
		for(int i=1;i<=n;i++)
		{
			for(int j=1;j<=n;j++)
			{
				a[i][j] = 0;
				for(int k=1;k<=n;k++)
				a[i][j] = (a[i][j]+b[i][k]*b[k][j]%mod)%mod;
			}
		}
	}
	void ksm()
	{
		while(m)
		{
			if(m&1) ksm1();
			ksm2();m >>= 1;
		}
	}
}mat;

signed main()
{
	IOS;
	cin>>n>>m;
	mat.a[1][1] = n-2,mat.a[1][2] = n-1;
	mat.a[2][1] = 1,mat.a[2][2] = 0;
	mat.n = 2;mat.m = m;mat.init();mat.ksm();
	va[1][1] = 0,va[2][1] = 1;
	for(int i=1;i<=2;i++)
	{
		for(int j=1;j<=1;j++)
		{
			for(int k=1;k<=2;k++)
			anw[i][j] = (anw[i][j]+mat.c[i][k]*va[k][j]%mod)%mod;
		}
	}
	cout<<anw[2][1];
	return 0;
}

可乐

题意：
就是给你n个城市，机器人在1号城市上，机器人有三种行为，停留在原地，去某个相邻的城市，或者自爆。机器人每秒都会行动一次，第0秒在1号点。现在问你经过了t秒，机器人的行为方案树是多少。

思考：

刚看到这题的时候没有太多的想法，就感觉这个转移的状态不少，然后感觉应该挺复杂的。但是一看时间1e6,点就30个。如果两重循环dp的话，3e7也还可以，不过空间卡的也很
转移的话那么相邻城市不用说，停留的话就自己想自己连边，自爆？其实可以设置一个自爆点0或者n+1,每个点连自爆点，这样就处理了自爆。那么最后的答案就是在t秒后每个点停留的方案数。不过当时有点疑惑的就是，对于自爆点是否需要自己连自己，因为自爆不能连续自爆，如果连边了会不会重复算，实际上不会，只是把这一个的时间点的自爆传到下一个时间点。当然也可以不自爆连边，把每个时间点自爆的点都加一次就可以了。
不过如果题目范围再大一点，这个dp就超时了。发现如果我把邻接矩阵写出来，每次转移是什么？就是自己乘以自己。如果转移t次呢？那么就是转移矩阵的t次方。那答案呢？答案就是转移后的1到每个点的方案数。

代码：

矩阵
#include<bits/stdc++.h>
#define fi first
#define se second
#define pb push_back
#define db double
#define int long long
#define PII pair<int,int >
#define mem(a,b) memset(a,b,sizeof(a))
#define IOS std::ios::sync_with_stdio(false),cin.tie(0),cout.tie(0);

using namespace std;
const int mod = 2017,inf = 1e18;
const int N = 2e5+10,M = 100;

int T,n,m,k;
int va[N];

struct matx{
	int n,m,a[M][M],b[M][M],c[M][M];
	void init()
	{
		for(int i=1;i<=n;i++)
		for(int j=1;j<=n;j++) c[i][j] = (i==j);
	}
	void ksm1()
	{
		for(int i=1;i<=n;i++)
		{
			for(int j=1;j<=n;j++)
			b[i][j] = c[i][j];
		}
		for(int i=1;i<=n;i++)
		{
			for(int j=1;j<=n;j++)
			{
				c[i][j] = 0;
				for(int k=1;k<=n;k++)
				c[i][j] = (c[i][j]+a[i][k]*b[k][j]%mod)%mod;
			}
		}
	}
	void ksm2()
	{
		for(int i=1;i<=n;i++)
		{
			for(int j=1;j<=n;j++)
			b[i][j] = a[i][j];
		}
		for(int i=1;i<=n;i++)
		{
			for(int j=1;j<=n;j++)
			{
				a[i][j] = 0;
				for(int k=1;k<=n;k++)
				a[i][j] = (a[i][j]+b[i][k]*b[k][j]%mod)%mod;
			}
		}
	}
	void ksm()
	{
		while(m)
		{
			if(m&1) ksm1();
			ksm2();m >>= 1;
		}
	}
}mat;

signed main()
{
	IOS;
	cin>>n>>m;
	for(int i=1;i<=m;i++)
	{
		int x,y;
		cin>>x>>y;
		mat.a[x][y] = mat.a[y][x] = 1;
	}
	for(int i=1;i<=n+1;i++) mat.a[i][i] = 1;
	for(int i=1;i<=n;i++) mat.a[i][n+1] = 1;
	cin>>k;
	mat.n = n+1,mat.m = k;mat.init();mat.ksm(); //注意是n+1个点，因为多了一个自爆点。
	int ans = 0;
	for(int i=1;i<=n+1;i++) ans = (ans+mat.c[1][i])%mod;
	cout<<ans;
	return 0;
}

直接dp
#include<bits/stdc++.h>
#define fi first
#define se second
#define pb push_back
#define db double
#define ll long long
#define PII pair<int,int >
#define mem(a,b) memset(a,b,sizeof(a))
#define IOS std::ios::sync_with_stdio(false),cin.tie(0),cout.tie(0);

using namespace std;
const int mod = 2017,inf = 1e18;
const int N = 1e6+1,M = 31;

int T,n,m,k;
int dp[N][M];

vector<int > e[M];

signed main()
{
	IOS;
	cin>>n>>m;
	for(int i=1;i<=m;i++)
	{
		int a,b;
		cin>>a>>b;
		e[a].pb(b);
		e[b].pb(a);
	}
	for(int i=1;i<=n;i++)
	{
		e[i].pb(0);
		e[i].pb(i);
	}
	e[0].pb(0);
	cin>>k;
	dp[0][1] = 1;
	for(int i=0;i<=k;i++)
	{
		for(int a=0;a<=n;a++)
		{
			for(auto b:e[a])
			dp[i+1][b] = (dp[i+1][b]+dp[i][a])%mod;
		}
	}
	int ans = 0;
	for(int i=0;i<=n;i++) ans = (ans+dp[k][i])%mod;
	cout<<ans;
	return 0;
}