高斯约旦消元法求逆矩阵的思想(分块矩阵的逆矩阵)

2022-07-29 13:20:56 浏览数 (1)

大家好,又见面了,我是你们的朋友全栈君。

luogu P4783 【模板】矩阵求逆

题目描述

求一个 N × N N×N N×N的矩阵的逆矩阵。答案对 1 0 9 7 10^9 7 109 7取模。

1.逆矩阵的定义

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

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

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

0.高斯-约旦消元

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

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

代码语言:javascript复制
高斯消元 -> 消成上三角矩阵 

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

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

代码语言:javascript复制
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("%.2lfn",a[i][n 1]/a[i][i]);
}

1.矩阵求逆

思路

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

原理 A − 1 ∗ [ A I ] = [ I A − 1 ] A^{-1} * [AI] = [I A^{-1}] A−1∗[AI]=[IA−1]

举个栗子 求 [ 2 − 1 0 − 1 2 − 1 0 − 1 2 ] left[ begin{matrix} 2 &amp; -1 &amp; 0 \ -1 &amp; 2 &amp; -1 \ 0 &amp; -1 &amp; 2 end{matrix} right] ⎣⎡​2−10​−12−1​0−12​⎦⎤​

首先 [ 2 − 1 0 1 0 0 − 1 2 − 1 0 1 0 0 − 1 2 0 0 1 ] begin{bmatrix} 2 &amp; -1 &amp; 0 &amp; 1 &amp; 0 &amp; 0 \ -1 &amp; 2 &amp; -1 &amp; 0 &amp; 1 &amp; 0 \ 0 &amp; -1 &amp; 2 &amp; 0 &amp; 0 &amp; 1 end{bmatrix} ⎣⎡​2−10​−12−1​0−12​100​010​001​⎦⎤​ 对左边进行消元可得 [ 2 − 1 0 1 0 0 0 3 2 − 1 1 2 1 0 0 0 4 3 1 3 2 3 1 ] left[ begin{matrix} 2 &amp; -1 &amp; 0 &amp; 1 &amp; 0 &amp; 0 \ 0 &amp; frac{3}{2} &amp; -1 &amp; frac{1}{2} &amp; 1 &amp; 0 \ 0 &amp; 0 &amp; frac{4}{3} &amp; frac{1}{3} &amp; frac{2}{3} &amp; 1 end{matrix} right] ⎣⎡​200​−123​0​0−134​​121​31​​0132​​001​⎦⎤​ 此时已消成上三角矩阵,高斯消元开始回代,但约旦会消成对角矩阵 [ 2 0 0 3 2 1 1 2 0 3 2 0 3 4 3 2 3 4 0 0 4 3 1 3 2 3 1 ] left[ begin{matrix} 2 &amp; 0 &amp; 0 &amp; frac{3}{2} &amp; 1 &amp; frac{1}{2} \ 0 &amp; frac{3}{2} &amp; 0 &amp; frac{3}{4} &amp; frac{3}{2} &amp; frac{3}{4} \ 0 &amp; 0 &amp; frac{4}{3} &amp; frac{1}{3} &amp; frac{2}{3} &amp; 1 end{matrix} right] ⎣⎡​200​023​0​0034​​23​43​31​​123​32​​21​43​1​⎦⎤​ 最后每行除以系数 [ 1 0 0 3 4 1 2 1 4 0 1 0 1 2 1 1 2 0 0 1 1 4 1 2 3 4 ] left[ begin{matrix} 1 &amp; 0 &amp; 0 &amp; frac{3}{4} &amp; frac{1}{2} &amp; frac{1}{4} \ 0 &amp; 1 &amp; 0 &amp; frac{1}{2} &amp; 1 &amp; frac{1}{2} \ 0 &amp; 0 &amp; 1 &amp; frac{1}{4} &amp; frac{1}{2} &amp; frac{3}{4} end{matrix} right] ⎣⎡​100​010​001​43​21​41​​21​121​​41​21​43​​⎦⎤​ 此时右半边即为所求

2.细节

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

C o d e Code Code

代码语言:javascript复制
#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("%lldn",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;
}

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

有什么问题欢迎评论区指出 :)

参考文章

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

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

发布者:全栈程序员栈长,转载请注明出处:https://javaforall.cn/129183.html原文链接:https://javaforall.cn

0 人点赞