P4783 【模板】矩阵求逆

传送门

线性代数真的好珂怕……以下如果有漏洞欢迎指出

定义矩阵的三种初等行变换:

1.交换某两行

2.将某一行的所有元素乘上\(k\)(\(k\neq 0\))

3.将某一行的所有元素乘上\(k\)加到另一行去

每一个初等变换都对应一个初等矩阵,即矩阵\(A\)做某一线性变换等价于用一个对应的初等矩阵左乘\(A\)。若有一堆初等变换\(1,2,...l\),对应的初等矩阵分别为\(P_1,P_2,...,P_l\),那么经过这些线性变换后的矩阵即为\(P_l....P_2P_1A=PA\)\(P\)为之前那堆东西的乘积)

对于一个矩阵\(A\)\(A\)可逆的充分必要条件是\(A\)经过若干次初等行变换可以变成\(E\)\(E\)即单位矩阵),即存在一个矩阵\(P\)使得\(PA=E\),则\(P=A^{-1}\)

通过初等行变换使得\(A\)变为\(E\)并不困难,可以用高斯消元解决,先消成上三角矩阵,然后再消成对角矩阵

考虑怎么求出\(P\),因为有\(PA=E,PE=P\),如果我们同时维护两个矩阵\(A,B\),令\(B\)一开始时等于\(E\),在把\(A\)变为\(E\)的过程中对\(B\)也做相等的初等变换,那么当\(A\)变为\(E\)时,\(B\)也就变为了\(P\)(因为做初等行变换等价于被对应的初等矩阵左乘)

如果在高斯消元的过程中发现无法将\(A\)变为\(E\),输出无解即可

//minamoto
#include<bits/stdc++.h>
using namespace std;
#define getc() (p1==p2&&(p2=(p1=buf)+fread(buf,1,1<<21,stdin),p1==p2)?EOF:*p1++)
char buf[1<<21],*p1=buf,*p2=buf;
int read(){
    int res,f=1;char ch;
    while((ch=getc())>'9'||ch<'0')(ch=='-')&&(f=-1);
    for(res=ch-'0';(ch=getc())>='0'&&ch<='9';res=res*10+ch-'0');
    return res*f;
}
char sr[1<<21],z[20];int C=-1,Z=0;
inline void Ot(){fwrite(sr,1,C+1,stdout),C=-1;}
void write(int x){
    if(C>1<<20)Ot();if(x<0)sr[++C]='-',x=-x;
    while(z[++Z]=x%10+48,x/=10);
    while(sr[++C]=z[Z],--Z);sr[++C]=' ';
}
const int N=405,mod=1e9+7;
int n;
struct Matrix{
    int a[N][N];
    inline void clr(){memset(a,0,sizeof(a));}
    int* operator [](int x){return a[x];}
    void SWAP(int x,int y){for(int i=1;i<=n;++i)swap(a[x][i],a[y][i]);}
    //交换某两行
    void MUL(int x,int k){for(int i=1;i<=n;++i)a[x][i]=(1ll*a[x][i]*k%mod+mod)%mod;}
    //将某一行的所有元素乘上k
    void MD(int x,int y,int k){for(int i=1;i<=n;++i)a[x][i]=((a[x][i]+(1ll*a[y][i]*k%mod))%mod+mod)%mod;}
    //将某一行的所有元素乘上k加到另一行去
	void print(){
        for(int i=1;i<=n;++i){
            for(int j=1;j<=n;++j)write(a[i][j]);
            sr[++C]='\n';
        }
    }
}A,B;
int ksm(int a,int b=(mod-2)){
    int res=1;
    for(;b;b>>=1,a=1ll*a*a%mod)if(b&1)res=1ll*res*a%mod;
    return res;
}
int main(){
//	freopen("testdata.in","r",stdin);
    n=read();
    for(int i=1;i<=n;++i)for(int j=1;j<=n;++j)A[i][j]=read();
    for(int i=1;i<=n;++i)B[i][i]=1;
    for(int i=1;i<=n;++i){
    	//消成上三角矩阵 
        if(!A[i][i]){
            for(int j=i+1;j<=n;++j)if(A[j][i]){
                A.SWAP(i,j),B.SWAP(i,j);break;
            }
        }
        if(!A[i][i])return puts("No Solution"),0;
        //如果消着消着某一列没有数了,说明无解 
        B.MUL(i,ksm(A[i][i])),A.MUL(i,ksm(A[i][i]));
        for(int j=i+1;j<=n;++j)
        B.MD(j,i,-A[j][i]),A.MD(j,i,-A[j][i]);
    }
    //消成对角矩阵 
    for(int i=n-1;i;--i)for(int j=i+1;j<=n;++j)
    B.MD(i,j,-A[i][j]),A.MD(i,j,-A[i][j]);
    B.print();return Ot(),0;
}
posted @ 2018-11-13 20:11  bztMinamoto  阅读(178)  评论(0编辑  收藏  举报
Live2D