矩阵求逆

其实这玩意去年也搞过不过就是TLE鹅已
我们知道如果\(ab=1\),则\(b\)\(a\)的逆元,那我们现在有两个矩阵\(A\),\(A^{-1}\),已知\(AA^{-1}=E\),则\(A^{-1}\)\(A\)的逆元
那么我们应该怎么求\(A{-1}\)呢?
如果我们用手算,那么可以先搞出来伴随矩阵,然后再用行列式除以\(A\)的行列式(这里建议百度百科)大概还会补上吧.....
当然还有另一种方法。
\(A\)乘它的逆矩阵得到单位矩阵,那么我们可以将\(A\)经过一系列操作把它变成单位矩阵。同时对另一个单位矩阵进行相同的操作,这样最终这个单位矩阵就会变成\(A\)的逆矩阵
为什么呢?
接下来内容请感性李姐理解
先康一康这个等式:\(AA^{-1}=E\)
我们每次对\(A\)进行操作,可看做是让\(A\)左乘或者右乘某个矩阵,同时等号右边的单位矩阵也乘相同的矩阵,当\(A\)变成单位矩阵的时候,设此时右边的矩阵是\(B\),则有\(EA^{-1}=B\),所以此时原来的单位矩阵就是现在的逆矩阵
我们可以通过高斯消元来把\(A\)变成单位矩阵
这里简单说一下高斯消元。
首先,我们先将矩阵消成一个上三角(当然下三角也是可以的)。然后从\((n-1,n)\)开始,\(A_{ik}\)-=\(A_{jk}\cdot A_{ik}\),其中\(j\)为当前要消去的元素所在的列,\(i\)表示当前元素所在的行。因为单位矩阵第\(j\)行第\(j\)列才是1.注意是一列一列的消
看看代码吧(蒟蒻放弃卡常于是乎吸了氧气)

#include<iostream>
#include<cstdio>
#include<cstring>
#include<cstdlib>
#include<algorithm>
#include<cmath>
#include<queue>
#include<vector>
#define pa pair<int,int>
using namespace std;
typedef long long ll;
ll read()
{
	char ch=getchar();
	ll x=0;bool f=0;
	while(ch<'0'||ch>'9')
	{
		if(ch=='-') f=1;
		ch=getchar();
	}
	while(ch>='0'&&ch<='9')
	{
		x=(x<<3)+(x<<1)+(ch^48);
		ch=getchar();
	}
	return f?-x:x;
} 
const int mod=1000000007;
int n;
struct jz{
	ll s[409][409];
}a,b;
void Swap(jz &nw,int st,int en)
{
	for(int i=1;i<=n;i++)
	 swap(nw.s[st][i],nw.s[en][i]);
}
ll ksm(ll a)//记得开longlong
{
	ll r=1;ll b=mod-2;
	while(b)
	{
		if(b&1) r=(r*a+mod)%mod;
		a=(a*a+mod)%mod;
		b>>=1;
	}
	return r;
}
void cheng(jz &nw,int ln,ll ny)
{
	for(int i=1;i<=n;i++)
	 nw.s[ln][i]=((nw.s[ln][i]*ny)%mod+mod)%mod;
}
void xiao(jz &nw,int nl,int fl,int v)
{
	for(int i=1;i<=n;i++)
	 nw.s[nl][i]=((nw.s[nl][i]+nw.s[fl][i]*v)%mod+mod)%mod;
}
void print()
{
		for(int i=1;i<=n;i++)
	{
		for(int j=1;j<=n;j++)
		 printf("%lld ",b.s[i][j]);
		printf("\n"); 
	}
	printf("\n");
}
int main()
{
	n=read();
	memset(a.s,0,sizeof(a.s));
	memset(b.s,0,sizeof(b.s));
	for(int i=1;i<=n;i++)
	 for(int j=1;j<=n;j++)
	  a.s[i][j]=read();
	for(int i=1;i<=n;i++)
	 for(int j=1;j<=n;j++)
	  b.s[i][i]=1;
	for(int i=1;i<=n;i++)//消成上三角矩阵
	{
		if(!a.s[i][i])
		 for(int j=i;j<=n;j++)
			if(a.s[j][i])
			{Swap(a,i,j),Swap(b,i,j);break;}	
		if(!a.s[i][i]){printf("No Solution\n");return 0;}//让对角线上不是0
	        ll ny=ksm(a.s[i][i]);
		cheng(a,i,ny);cheng(b,i,ny);//整行变换
		for(int j=i+1;j<=n;j++)//消去该列剩下的数
	        { xiao(b,j,i,-a.s[j][i]);
	          xiao(a,j,i,-a.s[j][i]);	
		}  
	}
	for(int i=n-1;i>=1;i--)//把上三角矩阵消成单位矩阵
	{
		for(int j=i+1;j<=n;j++)
		{   
		    xiao(b,i,j,-a.s[i][j]);
		    xiao(a,i,j,-a.s[i][j]);
		}
	}  
    print();
}
posted @ 2020-06-18 18:44  千载煜  阅读(448)  评论(0编辑  收藏  举报