高斯消元法

模板链接：https://www.luogu.org/problemnew/show/P3389

【模板】高斯消元法

题目背景

$\mathcal{Gauss}$消元

题目描述

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

输入输出格式

输入格式：

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

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

输出格式：

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

如果不存在唯一解，在第一行输出”No Solution”.

输入输出样例

输入样例#1：

3
1 3 4 5
1 4 7 3
9 3 2 2

输出样例#1：

-0.97
5.18
-2.39

说明

$1 \leq n \leq 100, \left | a_i \right| \leq {10}^4 , \left |b \right| \leq {10}^4$

题解

高斯消元法是用来解多元多次方程组的，我们将方程组的各个方程表示为矩阵形式，具体如下：

以下列方程组为例：

$\ X_1+X_2+X_3=0$
$2X_1+2X_2+X_3=1$
$\ X_1+X_2-X_3=-1$

写成矩阵形式如下：  

$\begin{matrix} 1&1&1&0\newline 2&2&1&1\newline 1&1&-1&-1\newline \end{matrix}$

之后，我们每次都考虑用一个方程去消元：

我们用第一个方程消掉第一个未知数，消元后的矩阵如下：

$\begin{matrix} 0&-1&1\newline   0&-2&-1\newline   \end{matrix}$

这样消元后，方程数-1，未知数数-1。

在消第二个未知数时，我们发现$x_2$的系数都为0，我们称这样的未知数为自由元。因为系数都为0，我们跳过对$x_2$的消元，矩阵变成下面这样：

$\begin{matrix} -1&1\newline   -2&-1\newline   \end{matrix}$

经过这样消元，未知数的数量-1。

最后一步，我们消去$x_3$：

$\begin{matrix} 0\newline   \end{matrix}$

矩阵只剩下0，这样这个方程就是有解的，反之无解。另外，如果有自由元的话这个方程就是有无穷组解的。当没有自由元时该方程有唯一解，此时我们可以将结果逐行代回求解。

代码实现方面，我们只需要模拟以上过程消元即可，单次消元的复杂度是$O(nm)$的，总消元次数为$O(min(n,m))$，所以总复杂度为$O(min(n,m)nm)$

代码

每次我们用系数最大的方程消元，利于减小误差。

#include<bits/stdc++.h>
#define eps 1e-8
using namespace std;
const int M=105;
double x[M][M];
int n;
void in()
{
    scanf("%d",&n);
    for(int i=0;i<n;++i)
    for(int j=0;j<=n;++j)
    scanf("%lf",&x[i][j]);
}
void ac()
{
    int pos;
    double ra;
    for(int i=0;i<n;++i)
    {
        pos=i;
        for(int j=i+1;j<n;++j)
        if(fabs(x[pos][i])<fabs(x[j][i]))pos=j;
        if(fabs(x[pos][i])<eps)
        {
            printf("No Solution\n");
            exit(0);
        }
        if(pos!=i)
        for(int j=0;j<=n;++j)swap(x[pos][j],x[i][j]);
        for(int j=i+1;j<n;++j)
        {
            ra=x[j][i]/x[i][i];
            for(int k=i;k<=n;++k)
            x[j][k]-=ra*x[i][k];
        }
    }
    for(int i=n-1;i>=0;--i)
    {
        for(int j=i+1;j<n;++j)x[i][n]-=x[j][n]*x[i][j];
        x[i][n]/=x[i][i];
    }
    for(int i=0;i<n;++i)
    printf("%.2lf\n",x[i][n]);
}
int main()
{
    in();ac();
    return 0;
}