高斯消元法

前置数学知识

n元线性方程是具有如下形式的方程：
$$a_1{x}_1+a_2{x}_2+a_3{x}_3+\dots +a_n{x}_n=b$$
其中，$a_1,a_2,...$以及常数项$b$均为已知的实数或者复数；$x_1,x_2,...x_n$是未知变量。$n$为正整数。
同理，一个$n$元线性方程组就是由多个$n$元线性方程组成的。一个含m个方程的n元线性方程组形如：
$$\{\begin{array}{c}{a}_{11}{x}_1+a_{12}{x}_2+...+a_{1n}{x}_n=b_1\\ a_{21}{x}_1+a_{22}{x}_2+...+a_{2n}{x}_n=b_2\\ ...\\ a_{m1}{x}_1+a_{m2}{x}_2+...+a_{mn}{x}_n=b_n\end{array}$$
其中，$a_{ij}$和$b_i$都是已知的数，我们称上诉方程组为$m\times n$线性方程组。

当$b_i$全为$0$时，称上诉线性方程组为齐次线性方程组。
当$b_i$不全为$0$时，称上诉线性方程组为非齐次线性方程组。

线性方程组的解有三种情况：（1）无解，（2）有唯一解，（3）有无穷多个解。如果一个线性方程组有唯一解或无穷多个解，则称它是相容的，如果无解则称它是不相容的。

高斯消元中的实现（线性代数部分）还涉及到矩阵相关的知识，这里也补充一下：
由$m\times n$个（实或复）数$a_{ij}(i=1,2,...,m;j=1,2,...,n)$排成$m$行$n$列的数表
$$\left[\begin{array}{cccc}{a}_{11}& a_{12}& \cdots & a_{1n}\\ a_{21}& a_{22}& \cdots & a_{2n}\\ ⋮& ⋮& \ddots & ⋮\\ a_{m1}& a_{m2}& \cdots & a_{mn}\end{array}\right]$$
称为$m$行$n$列矩阵（或$m\times n$矩阵）,简称矩阵。数$a_{ij}$称为矩阵的第$i$行第$j$列的元素。
一个线性方程组中每个变量的系数以及常数项按照它们在方程组中的位置构成的矩阵，称为方程组的增广矩阵。利用矩阵记号，我们可以简化求解线性方程组。

来个例子引入一下

$$\{\begin{array}{c}{x}_1-2x_2+x_3=0\\ 2x_2+x_3=7\\ 2x_1+x_2-x_3=1\end{array}$$
我们采用增广矩阵来记录对应的每次消元
原增广矩阵：
$$\left[\begin{array}{cccc}1& -2& 1& 0\\ 0& 2& 1& 7\\ 2& 1& -1& 1\end{array}\right]$$
保留第一个方程中的$x_1$，并且消去其余方程中的$x_1$
$$\left[\begin{array}{cccc}1& -2& 1& 0\\ 0& 2& 1& 7\\ 0& 5& -3& 1\end{array}\right]$$
保留第二个方程中的$x_2$，并且消去其余方程中的$x_2$
$$\left[\begin{array}{cccc}1& 0& 2& 7\\ 0& 2& 1& 7\\ 0& 0& -\frac{11}{2}& -\frac{33}{2}\end{array}\right]$$
保留第三个方程中的$x_3$，并且消去其余方程中的$x_3$，（同时将第三个方程的$x_3$的系数化为$1$）
$$\left[\begin{array}{cccc}1& 0& 0& 1\\ 0& 2& 0& 4\\ 0& 0& 1& 3\end{array}\right]$$
整理一下变成
$$\left[\begin{array}{cccc}1& 0& 0& 1\\ 0& 1& 0& 2\\ 0& 0& 1& 3\end{array}\right]$$
所以最后方程组的解为
$$\{\begin{array}{c}{x}_1=1\\ x_2=2\\ x_3=3\end{array}$$
其实，上诉所用的方法就叫高斯消元法，这里运用了初等行变换。

代码实现

有一点不一样的是，这个代码最后反代时用了类似递推的思想（从后往前计算答案）

#include<bits/stdc++.h>
#define re register
#define maxn 105
#define eps 1e-6
using namespace std;

int n;
double a[maxn][maxn],ans[maxn];

int main()
{
	ios::sync_with_stdio(false);
	cin.tie(0);
	cout.tie(0);
	cin >> n;
	for(re int i = 1;i <= n;i++)
		for(re int j = 1;j <= n+1;j++)
			cin >> a[i][j];
	for(re int i = 1;i <= n;i++)
	{
		int r = i;
		for(re int j = i+1;j <= n;j++)
			if(fabs(a[r][i]) < fabs(a[j][i]))r = j;
		if(fabs(a[r][i]) < eps)
		{
			cout << "No Solution" << '\n';
			return 0;
		}
		if(r != i)swap(a[r],a[i]);
		double div = a[i][i];
		for(re int j = i;j <= n+1;j++)
			a[i][j] = a[i][j]/div;
		for(re int j = i+1;j <= n;j++)
		{
			div = a[j][i];
			for(re int k = i;k <= n+1;k++)
				a[j][k] = a[j][k]-a[i][k]*div;
		}
	}
	ans[n] = a[n][n+1];
	for(re int i = n-1;i >= 1;i--)
	{
		ans[i] = a[i][n+1];
		for(re int j = i+1;j <= n;j++)
			ans[i] = ans[i]-a[i][j]*ans[j];
	}
	for(re int i = 1;i <= n;i++)
		cout << fixed << setprecision(2) << ans[i] << '\n';
	return 0;
}