[总结] 矩阵

矩阵乘法

今天学习了矩阵，~~针不辍~~。

定义（想不清楚回到定义）

矩阵乘法便是两个矩阵的指定行与列的每个元素分别相乘再求和。

定义 $ * $ 的代码：

mat operator *(const mat &x){
		mat ans;ans.init();
		for(int i=1;i<=8;i++)
			for(int j=1;j<=8;j++)
				for(int k=1;k<=8;k++)
					(ans.a[i][j]+=a[i][k]*x.a[k][j])%=P;
		return ans;			
}

Updated on 21th February 定义 * 的时候一定要写 $x.$。

Updated on 25th February 今天挨了模拟赛题的坑：

论矩阵的标准写法

struct mat{
	int a[maxn][maxn],n,m;
	inline void init(int _m = 0, int _n = 0) {
    	m = _m; n = _n;
	}
	mat(){
		memset(a,0,sizeof a);
	}
	mat operator *(const mat &x){
        mat ans;ans.init(m,x.n);
        for(int i=1;i<=m;i++)
            for(int j=1;j<=x.n;j++)
                for(int k=1;k<=n;k++)
                    (ans.a[i][j]+=1LL*a[i][k]*x.a[k][j])%=P;
        return ans;
    }
    mat power(int b){
    	mat ans;ans.init(n,n); 
    	for(int i=1;i<=n;i++)ans.a[i][i]=1;
    	mat tmp=(*this);
    	while(b){
    		if(b&1)ans=ans*tmp;
			tmp=tmp*tmp;
			b>>=1;
		}
		return ans;
	}
};

一定要把 $m$ 和 $n$ 记录，避免在有些题出现多余的转移。
快速幂写法，改掉递归写法的坏毛病，$tmp$ 要指向 $this$ 。

lrj都是什么年代的人物了
-gyx

平时写题不要省事，知道以后会出锅的地方平时就更应该注意。

以上为更新部分，白丢200pts 的经验教训

定义了乘号就可以直接将矩阵相乘了。

细节部分

毒瘤

每次一定要初始化矩阵元素都为 0 。

memset(a,0,sizeof a);

不然你会 wa 的很惨。

与图论的结合：

设图初始时的邻接矩阵为 $A_1$ ，那么矩阵是支持幂运算的。

$A_i$ 在矩阵上表示的是任意两点间的包含 i 条边的路径条数。

所以需要用到快速幂来表示$A^i$：

mat power(mat a,int b){
	if(!b)return base;//base是单位矩阵，主对脚线都是 1。
	mat k=power(a,b/2);
	k=k*k;
	if(b&1)k=k*a;
	return k;
}

例题：

例 1 P2233 公交车路线

和板子差不多，只不过让 $E$ 点的出边改为 0 ，因为到了 $E$ 就不能走了。

#include <iostream>
#include <cstdio>
#include <cstring>
using namespace std;
const int P = 1000;
struct mat{
	int a[9][9];
	mat init(int op=0){
		memset(a,0,sizeof a);
		if(op==1){
			for(int i=1;i<=8;i++)a[i][i]=1;
		}
		if(op==2){
			for(int i=1;i<=7;i++)a[i][i+1]=a[i+1][i]=1;
			a[5][4]=a[5][6]=0;
			a[8][1]=a[1][8]=1;
		}
	}
	mat operator *(const mat &x){
		mat ans;ans.init();
		for(int i=1;i<=8;i++)
			for(int j=1;j<=8;j++)
				for(int k=1;k<=8;k++)
					(ans.a[i][j]+=a[i][k]*x.a[k][j])%=P;
		return ans;			
	}
}A,base;
mat power(mat a,int b){
	if(!b)return base;
	mat k=power(a,b/2);
	k=k*k;
	if(b&1)k=k*a;
	return k;
}
int n;
int main(){
	A.init(2);base.init(1);
	scanf("%d",&n);
	printf("%d",power(A,n).a[1][5]);
	return 0;
}

Updated on 21th,February

今天终于把线性代数 A 完了，来补一补例题。

例 2 . P1397 NOI2013 矩阵游戏

NOI签到题qwq

题意分析

给你一个 $n * m $ 的矩阵和两种构造方式:

$F[1][1]=1$

$F[i,j]=a * F[i,j-1]+b$

$F[i,1]=c * F[i-1,m]+d$

题解

无非是序列的两种递推方式，因此可以找到序列的转移矩阵。

$A_1$：

$\left(\begin{matrix}f_{i,j} & 1 \end{matrix}\right) = \left(\begin{matrix} f_{i,j-1} & 1 \end{matrix}\right)$ $\left(\begin{matrix}a&0\\b&1\end{matrix}\right)$

$A_2$:

$\left(\begin{matrix} f_{i,1}& 1 \end{matrix}\right)=$
$\left(\begin{matrix}f_{i-1,m} & 1 \end{matrix}\right)$
$\left(\begin{matrix}c & 0\\d & 1 \end{matrix}\right)$

因此答案为 $A_1^{m-1}(A_2A_1^{m-1})^{n-1}T$

T是初始矩阵:

$\left(\begin{matrix}1 & 1\\0 & 0 \end{matrix}\right)$

另外需要高精度存储在转化为数字，可以边转化边取模。

细节

初始矩阵为了乘方便可以设成 2 * 2的，这并不影响结果。
取模的时候是根据矩乘的费马小定理，如果$a\ne1$，那么取模的时候相当于对 $P-1$ 取模。

#include <iostream>
#include <cstdio>
#include <cstring>
#define ll long long
using namespace std;
const int P = 1000000007;
ll a,b,c,d,n,m;
char s1[1000005],s2[1000005];
struct mat{
	ll a[3][3];
	mat(){
		memset(a,0,sizeof a);
	}
	mat operator *(const mat &x){
		mat ans;
		for(int i=1;i<=2;i++)
			for(int j=1;j<=2;j++)
				for(int k=1;k<=2;k++)
					(ans.a[i][j]+=a[i][k]*x.a[k][j])%=P;
	return ans;
	}
	
}A,B,T,ans,res;
void power(mat a,int b){
	ans.a[1][1]=ans.a[2][2]=1;
	ans.a[1][2]=ans.a[2][1]=0;
	while(b){
		if(b&1)ans=a*ans;
		a=a*a;
		b>>=1;
	}
}
int main(){
	scanf("%s%s%lld%lld%lld%lld",s1,s2,&a,&b,&c,&d);
	A.a[1][1]=a;A.a[1][2]=0;A.a[2][1]=b;A.a[2][2]=1;
	B.a[1][1]=c;B.a[1][2]=0;B.a[2][1]=d;B.a[2][2]=1;
	T.a[1][1]=1;T.a[1][2]=1;
	ll p=P-1;
	if(a==1)p=P;
	int len=strlen(s1);
	for(int i=0;i<len;i++)n=(n*10+s1[i]-'0')%p;
	len=strlen(s2);
	for(int i=0;i<len;i++)m=(m*10+s2[i]-'0')%p;
	power(A,m-1);
	for(int i=1;i<=2;i++)
		for(int j=1;j<=2;j++)
			A.a[i][j]=res.a[i][j]=ans.a[i][j];
	A=A*B;
	power(A,n-1);
	res=ans*res;
	res=T*res;
	printf("%lld",res.a[1][1]);
	return 0;		
}

例 3 .P4159 SCOI2009 迷路

题意分析

如果这道题没有边权，那么就是一道板子题。

注意，邻接矩阵解决的题边权都应该是 1 的。

题解

考虑如何处理边权。

借用分层图的思想，对于一条路径 $(u,v,w)$ ，$u_0$是真正的点，其他下标的点都是辅助点。所以就可以在 $u_0$ 和 $v_{w-1}$ 之间连上一条权值为 1 的边。

可以这么理解：

你花 1 的长度走到了假的 $v$ 点，还需要向上一层一层爬，爬 $w-1$ 个距离爬到真正的 $v$ 点。

所以我们把每一层的对应点间连上一条权值为 1 的边。

这道题用到了等价转换思想，不管是whk还是竞赛中不可缺少的思想。

#include <iostream>
#include <cstdio>
#include <cstring>
#include <vector>
using namespace std;
const int P = 45989;
const int maxn = 150;
int n,m,T,s,t;
struct node{
	int to,nxt;
}e[maxn<<1];
int head[maxn],cnt=1;
void link(int u,int v){
	e[++cnt].to=v;e[cnt].nxt=head[u];head[u]=cnt;
}
struct mat{
	int a[maxn][maxn];
	mat(){
		memset(a,0,sizeof a);
	}
	mat operator *(const mat &x){
		mat ans;
		for(int i=1;i<=m;i++)
			for(int j=1;j<=m;j++)
				for(int k=1;k<=m;k++)
					(ans.a[i][j]+=a[i][k]*x.a[k][j])%=P;
		return ans;			
	}
}num,base,based;
mat power(mat a,int b){
	if(!b)return based;
	mat k=power(a,b/2);
	k=k*k;
	if(b&1)k=k*a;
	return k;
}
int main(){
	scanf("%d%d%d%d%d",&n,&m,&T,&s,&t);
	s++;t++;
	for(int i=1,u,v;i<=m;i++){
		scanf("%d%d",&u,&v);
		u++,v++;
		link(u,v);link(v,u);
	}
	for(int i=2;i<=cnt;i++){
		int u=e[i].to;
		for(int j=head[u];j;j=e[j].nxt){
			if(j!=(i^1))num.a[i][j]++;
		}
	}
	for(int i=head[s];i;i=e[i].nxt)base.a[1][i]++;
	for(int i=1;i<=cnt;i++)based.a[i][i]=1;
	m=cnt;
	num=base*power(num,T-1);
	/*for(int i=1;i<=m;i++)
		for(int j=1;j<=m;j++)cerr<<num.a[i][j]<<endl;*/
	long long ans=0;
	for(int i=head[t];i;i=e[i].nxt){
		ans=(ans+num.a[1][i^1])%P;
	}
	cout<<ans;
	return 0;
}

例 4 . P2151 SDOI2009 HH去散步

题意分析

如果可以回到原来的位置，那么这又是一个板子题。

题解

仔细想想好像没有什么好的处理方法。

这里用到一个非常重要的思想：点边互换。

什么意思呢，把图中的边作为新图的点建图，使之可以形成一个 $DAG$ 。

可以这么理解：（仅仅只是理解）

一般的邻接矩阵是从点开始按照距离扩展的。
所以如果把边看成点的话， $t=0.5$ 是初始位置， $t=k-0.5$ 是最终位置。

所以转化就转化完了，关键这还是一个代码实现难度比较高的题。

具体做法：

枚举每一条边并向外扩展，如果不是反边假邻接矩阵++。
可以有一个虚拟边 1 ( 因为cnt是从1开始，边的编号最小是二 )，初始一下距离为一可以到达的边集便是源点 s 的所有出边。
对矩阵进行快速幂操作，幂是 $T-1$ ，也就是初始状态定义好后进行 $T-2$次转移。

确实是一道特别有思维含量的题。

#include <iostream>
#include <cstdio>
#include <cstring>
#include <vector>
using namespace std;
const int P = 45989;
const int maxn = 150;
int n,m,T,s,t;
struct node{
	int to,nxt;
}e[maxn<<1];
int head[maxn],cnt=1;
void link(int u,int v){
	e[++cnt].to=v;e[cnt].nxt=head[u];head[u]=cnt;
}
struct mat{
	int a[maxn][maxn];
	mat(){
		memset(a,0,sizeof a);
	}
	mat operator *(const mat &x){
		mat ans;
		for(int i=1;i<=m;i++)
			for(int j=1;j<=m;j++)
				for(int k=1;k<=m;k++)
					(ans.a[i][j]+=a[i][k]*x.a[k][j])%=P;
		return ans;			
	}
}num,base,based;
mat power(mat a,int b){
	if(!b)return based;
	mat k=power(a,b/2);
	k=k*k;
	if(b&1)k=k*a;
	return k;
}
int main(){
	scanf("%d%d%d%d%d",&n,&m,&T,&s,&t);
	s++;t++;
	for(int i=1,u,v;i<=m;i++){
		scanf("%d%d",&u,&v);
		u++,v++;
		link(u,v);link(v,u);
	}
	for(int i=2;i<=cnt;i++){
		int u=e[i].to;
		for(int j=head[u];j;j=e[j].nxt){
			if(j!=(i^1))num.a[i][j]++;
		}
	}
	for(int i=head[s];i;i=e[i].nxt)base.a[1][i]++;
	for(int i=1;i<=cnt;i++)based.a[i][i]=1;
	m=cnt;
	num=base*power(num,T-1);
	/*for(int i=1;i<=m;i++)
		for(int j=1;j<=m;j++)cerr<<num.a[i][j]<<endl;*/
	long long ans=0;
	for(int i=head[t];i;i=e[i].nxt){
		ans=(ans+num.a[1][i^1])%P;
	}
	cout<<ans;
	return 0;
}

Updated on 2021/8/6

~~快退役了才发现还有个细节~~

矩阵快速幂加速递推应用细节。

如果说 $j$ 状态可以转移到 $i$ 状态，也就是 $j->i$。

如果用行向量，是正向的，即转移矩阵 $a[j][i]=1$。

如果用列向量，是反向的，即转移矩阵 $a[i][j]=1$。

posted @ 2021-08-12 16:40 ¶凉笙阅读(65) 评论(0) 编辑收藏举报

会员力量，点亮园子希望

刷新页面返回顶部

¶凉笙

[总结] 矩阵

矩阵乘法

定义（想不清楚回到定义）

论矩阵的标准写法

以上为更新部分，白丢200pts 的经验教训

细节部分

与图论的结合：

例题：

例 1 P2233 公交车路线

Updated on 21th,February

例 2 . P1397 NOI2013 矩阵游戏

题意分析

题解

细节

例 3 .P4159 SCOI2009 迷路

题意分析

题解

例 4 . P2151 SDOI2009 HH去散步

题意分析

题解

Updated on 2021/8/6

公告