矩阵求逆

· · 个人记录 

#include<cstdio>
#include<iostream>
#define re register
#define ll long long
#define maxn 505
using namespace std;
const int mod=1e9+7;
ll map[maxn][maxn<<1];
int n,m;
ll quickpow(ll a,int b)
{
    ll ans=1;
    while(b)
    {
        if(b&1)
        ans=(ans*a)%mod;
        a=(a*a)%mod;
        b>>=1;
    }
    return ans;
}
int main()
{
    scanf("%d",&n);
    for(re int i=1;i<=n;++i)
    {
        for(re int j=1;j<=n;++j)
        {
            scanf("%d",&map[i][j]);
            map[i][i+n]=1;
        }
    }
    m=2*n;
    for(re int i=1;i<=n;++i)
    {
        for(re int j=i;j<=n;++j)
        {
            if(map[j][i]) 
            {
                for(re int k=1;k<=m;++k)
                swap(map[j][k],map[i][k]);//跟精度无关，找到直接换 
                break;
            }
        }
        if(!map[i][i]){printf("No Solution\n");return 0;}
        ll r=quickpow(map[i][i],mod-2);//找到逆元 
        for(re int j=1;j<=m;++j) 
        map[i][j]=(map[i][j]*r)%mod;//整整一行都要乘，相当于单位矩阵 
        for(re int j=1;j<=n;++j)//最终要消成单位矩阵所以每一行都要消除，高斯消元上三角就行 
        {
            if(i==j) continue;//注意 
            r=map[j][i];
            for(re int k=i;k<=m;++k)
                map[j][k]=(map[j][k]-r*map[i][k]%mod+mod)%mod;
        }
    }
    for(re int i=1;i<=n;++i)
    {
        for(re int j=n+1;j<=m;++j)
        {
            printf("%d",map[i][j]);
            if(j==m) printf("\n");
            else printf(" ");
        }
    }
    return 0;
}