矩阵求逆 —— 初等变换法（高斯-约旦消元）

luogu P4783 【模板】矩阵求逆

题目描述

求一个$N\times N$的矩阵的逆矩阵。答案对$10^9+7$取模。

1.逆矩阵的定义

假设 $A$ 是一个方阵，如果存在一个矩阵 $A^{-1}$，使得
$$A^{-1}A=I$$
并且
$$A{A}^{-1}=I$$

那么，矩阵 A 就是可逆的，$A^{-1}$ 称为 A 的逆矩阵

2.逆矩阵求法 —— 初等变换法（高斯-约旦消元）

0.高斯-约旦消元

详见P3389 【模板】高斯消元法题解部分

高斯约旦消元与高斯消元区别：

高斯消元 -> 消成上三角矩阵 

高斯-约旦消元 -> 消成对角矩阵 

约旦消元法的精度更好,代码更简单,没有回带的过程

void Gauss_jordan(){
	/***** 行的交换&加减消元 *****/ 
	for(re int i=1,r;i<=n;++i){	//正在处理第i行 
		r=i;
		for(re int j=i+1;j<=n;++j) 
			if(fabs(a[j][i])>fabs(a[r][i])) r=j;
		if(fabs(a[r][i])<eps){
			puts("No Solution");return;
		}
		if(i!=r) swap(a[i],a[r]);
		
		for(re int k=1;k<=n;++k){
		//每一行都处理  
			if(k==i) continue;
			double p=a[k][i]/a[i][i];
			for(re int j=i;j<=n+1;++j) a[k][j]-=p*a[i][j];
		} 
	}	
	
	//上述操作后会剩下对角矩阵,答案要除以系数    
	for(re int i=1;i<=n;++i) printf("%.2lf\n",a[i][n+1]/a[i][i]);
}

1.矩阵求逆

思路

• 求$A$的逆矩阵，把$A$和单位矩阵$I$放在一个矩阵里
• 对$A$进行加减消元使$A$化成单位矩阵
• 此时原来单位矩阵转化成逆矩阵

原理
$$A^{-1}\ast [AI]=[I{A}^{-1}]$$

举个栗子
求
$$\left[\begin{array}{ccc}2& -1& 0\\ -1& 2& -1\\ 0& -1& 2\end{array}\right]$$

首先
$$\left[\begin{array}{cccccc}2& -1& 0& 1& 0& 0\\ -1& 2& -1& 0& 1& 0\\ 0& -1& 2& 0& 0& 1\end{array}\right]$$
对左边进行消元可得
$$\left[\begin{array}{cccccc}2& -1& 0& 1& 0& 0\\ 0& \frac{3}{2}& -1& \frac{1}{2}& 1& 0\\ 0& 0& \frac{4}{3}& \frac{1}{3}& \frac{2}{3}& 1\end{array}\right]$$
此时已消成上三角矩阵，高斯消元开始回代，但约旦会消成对角矩阵
$$\left[\begin{array}{cccccc}2& 0& 0& \frac{3}{2}& 1& \frac{1}{2}\\ 0& \frac{3}{2}& 0& \frac{3}{4}& \frac{3}{2}& \frac{3}{4}\\ 0& 0& \frac{4}{3}& \frac{1}{3}& \frac{2}{3}& 1\end{array}\right]$$
最后每行除以系数
$$\left[\begin{array}{cccccc}1& 0& 0& \frac{3}{4}& \frac{1}{2}& \frac{1}{4}\\ 0& 1& 0& \frac{1}{2}& 1& \frac{1}{2}\\ 0& 0& 1& \frac{1}{4}& \frac{1}{2}& \frac{3}{4}\end{array}\right]$$
此时右半边即为所求

2.细节

1. 开long long（不要冒风险，乘法很容易溢出）
2. 模意义下除以一个数等于乘上逆元，可用快速幂求逆元（费马小定理）

$Code$

#include<iostream>
#include<cstdio>
#include<cmath>
#define re register
#define il inline
#define ll long long
using namespace std;

il ll read(){
    ll s=0,f=0;char c=getchar();
    while(c<'0'||c>'9') f=(c=='-'),c=getchar();
    while(c>='0'&&c<='9') s=(s<<3)+(s<<1)+(c^'0'),c=getchar();
    return f?-s:s;
}

const int N=405,mod=1e9+7;
int n;
ll a[N][N<<1];
il ll qpow(ll x,ll k){
	ll ans=1;
	while(k){
		if(k&1) ans=ans*x%mod;
		x=x*x%mod;
		k>>=1;
	}
	return ans%mod;
}

il void Gauss_j(){	
	for(re int i=1,r;i<=n;++i){
		r=i;
		for(re int j=i+1;j<=n;++j)
			if(a[j][i]>a[r][i]) r=j;
		if(r!=i) swap(a[i],a[r]);
		if(!a[i][i]){puts("No Solution");return;}
		
		int kk=qpow(a[i][i],mod-2);	//求逆元 
		for(re int k=1;k<=n;++k){
			if(k==i) continue;
			int p=a[k][i]*kk%mod;
			for(re int j=i;j<=(n<<1);++j) 
				a[k][j]=((a[k][j]-p*a[i][j])%mod+mod)%mod;
		} 
		
		for(re int j=1;j<=(n<<1);++j) a[i][j]=(a[i][j]*kk%mod);
		//更新当前行 如果放在最后要再求一次逆元,不如直接放在这里  
	}	
	
	for(re int i=1;i<=n;++i){
		for(re int j=n+1;j<(n<<1);++j) printf("%lld ",a[i][j]);
		printf("%lld\n",a[i][n<<1]);
	}
}
int main(){
	n=read();
	for(re int i=1;i<=n;++i)
		for(re int j=1;j<=n;++j)
			a[i][j]=read(),a[i][i+n]=1;
	
	Gauss_j();
    return 0;
}

网上浏览一圈头都要炸掉，线性代数太可怕了，定义好多
最后只看懂了这种方法

有什么问题欢迎评论区指出 ：）

参考文章

线性代数之——矩阵乘法和逆矩阵

逆矩阵的几种求法与解析(很全很经典)