高斯消元 模板

传送门

给定一个n元一次方程组求解。

首先我们思考一下，自己解普通的二元一次方程组的时候是怎么解的？

我们肯定是要先用一个变量表示另一个变量，换句话说，先进行消元，然后转变为一元一次方程求解。

我们一般有两种做法，加减消元或者带入消元。（这里以二元一次方程组为例，多元的也一样）

但是我们在使用程序计算的时候需要一种有规律可寻的算法，所以我们选择在开始的时候进行加减消元，返回的时候进行带入消元（把计算出来的结果带回去再计算）

首先，我们获得了一个n元一次方程组。之后我们要消除一个元，我们会怎么做？先把这个式子的第一个元的系数变为1，再把这个式子所有未知数前面的系数都乘以一个实数，使得在和其他式子计算的时候把这个元的系数消为0。（脑补一下做题时候的画面）之后再诸如此种操作去不断消元，直到最后获得一个一元一次方程，就可以轻松的计算出最后的得数。

之后我们把这个得数回带到刚才的式子里面，就又可以解出一个得数（因为系数为1嘛），之后这样不断地返回就可以解出n元一次方程组的解。时间复杂度O（n^3）。

注意因为我们只能使用double，所以必然会有精度误差（除非强大到全写高精度……那是不可能的），所以我们尽量每次用系数最大的一个去消其他的元，这样能尽量减小误差。我们来偷一位dalao的证明：

这样高斯消元就结束了……（如果中间有不理解的直接想想自己平时解方程组怎么解就好了）

看一下代码。

#include<iostream>
#include<cstdio>
#include<cmath>
#include<algorithm>
#include<queue>
#include<cstring>
#include<utility>
#include<map>
#define pr pair<int,int>
#define mp make_pair
#define fi first
#define sc second
#define rep(i,a,n) for(int i = a;i <= n;i++)
#define per(i,n,a) for(int i = n;i >= a;i--)
#define enter putchar('\n')
using namespace std;
typedef long long ll;
const int M = 100005;
const int N = 500005;
const int INF = 1000000009;
const double eps = 1e-7;
int read()
{
    int ans = 0,op = 1;
    char ch = getchar();
    while(ch < '0' || ch > '9')
    {
    if(ch == '-') op = -1;
    ch = getchar();
    }
    while(ch >='0' && ch <= '9')
    {
    ans *= 10;
    ans += ch - '0';
    ch = getchar();
    }
    return ans * op;
}

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

int main()
{
    n = read();
    rep(i,1,n)
    rep(j,1,n+1) scanf("%lf",&g[i][j]);
    rep(i,1,n)
    {
    int r = i;
    rep(j,i+1,n) if(fabs(g[r][i]) < fabs(g[j][i])) r = j;//找到当前系数最大的
    if(fabs(g[r][i]) < eps)//如果系数为0说明无解
    {
        printf("No Solution");
        return 0;
    }
    if(i != r) swap(g[i],g[r]);//把两行交换，这样就相当于我们每次都消去第一行，比较方便
    double div = g[i][i];//取这个系数
    rep(j,i,n+1) g[i][j] /= div;//把这行的每个系数都/div（相当于把当前要消的元的系数变为1）
    rep(j,i+1,n)
    {
        div = g[j][i];
        rep(k,i,n+1) g[j][k] -= g[i][k] * div;//对每一个式子进行消元（其中i肯定被消了，因为原来的式子里其系数为1）
    }
    }
    ans[n] = g[n][n+1];//得到最终结果
    per(i,n-1,1)//开始回带
    {
    ans[i] = g[i][n+1];//取上一次的答案
    rep(j,i+1,n) ans[i] -= g[i][j] * ans[j];//计算新的答案，就是这一行的值-前面已经算出来的变量×系数，不理解可以自己模拟一下
    }
    rep(i,1,n) printf("%.2lf\n",ans[i]);//输出每个答案
    return 0;
}