• 矩阵求逆


    其实这玩意去年也搞过不过就是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("
    "); 
    	}
    	printf("
    ");
    }
    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
    ");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();
    }
    
  • 相关阅读:
    MongoDB
    Mac下将Python2.7改为Python3
    Kafka
    Server 基础概念 备忘
    应用内支付
    Sign in with apple
    Linux三剑客grep/sed/awk
    React-Native中使用到的一些JS特性
    Date Picker控件:
    Picker View控件:
  • 原文地址:https://www.cnblogs.com/lcez56jsy/p/13159192.html
Copyright © 2020-2023  润新知