高斯消元详解

题目描述

题目背景
Gauss消元

题目描述

给定一个线性方程组，对其求解

输入格式

第一行，一个正整数 $n$

第二至 $n+1$行，每行 $n+1$ 个整数，为$a_1,a_2\cdots a_n$和$b$，代表一组方程。

输出格式

共$n$行，每行一个数，第 $i$ 行为 $x_i$ （保留2位小数）

如果不存在唯一解，在第一行输出"$NoSolution$".

原题链接：P3389 【模板】高斯消元法

解题思路

我们先来研究研究这个方程组：
$$\{\begin{array}{l}{x}_1+2x_2-x_3=-6\\ 2x_1+x_2-3x_3=-9\\ -x_1-x_2+2x_3=7\end{array}$$
我们可以只看系数，构成一个系数矩阵，然后再加上等号右侧的常数，构成一个增广矩阵：

$\left(\begin{array}{cccc}1& 2& -1& -6\\ 2& 1& -3& -9\\ -1& -1& 2& 7\end{array}\right)$

用高斯消元解方程组的方法其实可以概括为增广矩阵的三类操作：

1. 用一个非0的数乘某一行。
2. 把其中一行的若干倍加到另一行上。
3. 交换两行的位置。

这三类操作便称为矩阵的初等行变换。

于是，高斯消元的大致思路就是：

1. 通过“矩阵的初等行变换”，消去系数，使得矩阵呈现三角形，如下：
   $\left(\begin{array}{cccc}1& 2& -1& -6\\ 0& 1& 1& 1\\ 0& 0& 1& 3\end{array}\right)$
   即第$n$行的前$n-1$个数都是$0$。
   这样的矩阵就叫做阶梯型矩阵。
   于是，我们可以解得$x_3=3$。
2. 由于第一步中得到了$x_3$的解，所以我们将这个解回代到上一层，便可以求出$x_2$。之后便是再将$x_3,x_2$回代，又可以求出$x_1$，于是矩阵可以化为：
   $\left(\begin{array}{cccc}1& 0& 0& 1\\ 0& 1& 0& -2\\ 0& 0& 1& 3\end{array}\right)$
   该矩阵给出了方程的解，叫做简化阶梯矩阵。

回代的过程比较简单，所以下面说一说如何进行消元。

来个具体的例子看看：
还是上文的方程组，我们对其进行如下操作：

• 我们通过一重循环，一列一列地将第$i$列的$i-1$个数消去为0。我们首先选择该列绝对值最大的数，并把该行换为第$i$行。其次，通过这行消去第$i+1$到$n$行，即将该列系数化为同一个数，再相减。如下：
  $\left(\begin{array}{cccc}1& 2& -1& -6\\ 2& 1& -3& -9\\ -1& -1& 2& 7\end{array}\right)$
  第一列绝对值最大的是2，将这一行换为第一行：
  $\left(\begin{array}{cccc}2& 1& -3& -9\\ 1& 2& -1& -6\\ -1& -1& 2& 7\end{array}\right)$
  与第2行相消，即将第1行的所有数乘上$\frac{1}{2}$,再减去第2行：
  $\left(\begin{array}{cccc}2& 1& -3& -9\\ 2\cdot \frac{1}{2}-1& 1\cdot \frac{1}{2}-2& -3\cdot \frac{1}{2}-(-1)& -9\cdot \frac{1}{2}-(-6)\\ -1& -1& 2& 7\end{array}\right)$
  同理，其他行也可以这么消去。
  接下来，看代码——

参考代码

#include<stdio.h>
#include<iostream>
#define eps 1e-10
using namespace std;

int n;
double g[105][105],ans[105];

double abs(double a) //手写好像更快
{
	return a>0?a:-a;
}

void Gauss()
{
	for(int c=1;c<=n;c++) //一列一列消元
	{
		int r=c;
		for(int i=c+1;i<=n;i++) //找绝对值最大数
			r=abs(g[i][c])>abs(g[r][c])?i:r;
		if(abs(g[r][c])<eps) //坑点！如果最大的数为0（或是接近0），说明该行已经都是0，也说明该方程无解。由于是double运算，可能结果很小但不是0，所以要设置精度
		{
			printf("No Solution\n");
			return;
		}
		if(r!=c) swap(g[c],g[r]); //换行
		double tmp=g[c][c];
		for(int j=c;j<=n+1;j++) g[c][j]/=tmp; //将该行该列系数化为1，方便计算
		for(int i=c+1;i<=n;i++) //相减，消元
		{
			tmp=g[i][c];
			for(int j=c;j<=n+1;j++) g[i][j]-=(g[c][j]*tmp);
		}
	}
	ans[n]=g[n][n+1];
	for(int i=n-1;i>=1;i--) //回代，输出答案
	{
		for(int j=n;j>i;j--) g[i][n+1]-=(ans[j]*g[i][j]);
		ans[i]=g[i][n+1];
	}
	for(int i=1;i<=n;i++) printf("%.2lf\n",ans[i]);
}

int main()
{
	scanf("%d",&n);
	for(int i=1;i<=n;i++)
		for(int j=1;j<=n+1;j++) scanf("%lf",&g[i][j]);
	Gauss();
	return 0;
}