约旦-高斯消元

约旦-高斯消元

本蒟蒻感觉这个容易理解一点，而且精度好像更优秀，还不用回代，就选择了这个

高斯消元的作用：把一个方程组化成上三角(高斯消元)或对三角(约旦-高斯消元)，使其容易求解

这里一化为对三角为例，其步骤如下

1. 选择一个没有被选过的未知数作为主元，选择一个包含这个主元的方程
2. 把这个主元系数化为 1
3. 利用即将学的加减消元，消去这个方程中的其他系数
4. 重复以上步骤

无解判定：消元时会用到除法，若除数为 0，则无解

求解：直接用唯一剩下的系数求就可以了

例：一个模板

#include<bits/stdc++.h>
using namespace std;
const double eps=1e-7;
int n,m;
double x[105][105];
int main() {
	scanf("%d",&n),m=n+1;
	for(int i=1;i<=n;i++)
		for(int j=1;j<=m;j++)
			scanf("%lf",&x[i][j]);
	for(int i=1,r;i<=n;i++) {
        // 找到最大的系数
		r=i;
		for(int j=i+1;j<=n;j++)
			if(fabs(x[r][i]-x[j][i])<=eps)
				r=j;
        // 无解判定
		if(fabs(x[i][i])<=eps)
			return printf("No Solution"),0;
		if(r^i)swap(x[i],x[r]);
        // 加减消元
		double tmp;
		for(int j=1;j<=n;j++)
			if(i^j) {
				tmp=x[j][i]/x[i][i];
				for(int k=i+1;k<=m;k++)
					x[j][k]-=x[i][k]*tmp;
			}
	}
	for(int i=1;i<=n;i++)
	x[i][m]/=x[i][i];
	for(int i=1;i<=n;i++)
		printf("%.2lf\n",x[i][m]);
}

逆矩阵

定义：假设 \(A\) 是一个方阵，如果存在 \(A^{-1}\) 使得 \(A^{-1}A=I\) 那么矩阵 \(A\) 可逆，\(A^{-1}\) 称为 \(A\) 的逆矩阵

思路：

1. 把 \(I\) 放在 \(A\) 的右边
2. 对 \(A\) 进行消元使得 \(A\) 称为单位矩阵，可直接使用约旦-高斯消元
3. 原来的单位矩阵转换成逆矩阵

因为：\(A^{-1}*[AI]=[IA^{-1}]\)

例：矩阵求逆。注意求逆元即可。

#include<bits/stdc++.h>
using namespace std;
typedef long long LL;
const LL mod=1000000007;
inline LL Pow(LL x,LL y) {
	register LL res=1;
	for(;y;y>>=1,x=x*x%mod)
	if(y&1)res=res*x%mod;
	return res;
}
int n,m;
LL x[405][805];
int main() {
	scanf("%d",&n),m=n<<1;
	for(int i=1;i<=n;i++)
		for(int j=1;j<=n;j++)
			scanf("%lld",&x[i][j]);
	for(int i=1;i<=n;i++)
		for(int j=1;j<=n;j++)
			x[i][j+n]=(i==j);
	for(int i=1,r;i<=n;i++) {
		r=i;
		for(int j=i+1;j<=n;j++)
			if(x[j][i]>x[r][i])r=j;
		if(!x[i][i])return puts("No Solution"),0;
		if(r^i)swap(x[i],x[r]);
		register LL inv=Pow(x[i][i],mod-2),tmp;
		for(int j=1;j<=n;j++)
			if(i^j) {
				tmp=x[j][i]*inv%mod;
				for(int k=i+1;k<=m;k++)
					x[j][k]=(x[j][k]-(x[i][k]*tmp%mod)+mod)%mod;
			}
		for(int j=1;j<=m;j++)
			x[i][j]=x[i][j]*inv%mod;
	}
	for(int i=1;i<=n;i++,putchar(10))
		for(int j=1;j<=n;j++)printf("%lld ",x[i][j+n]);
}