高斯消元

Posted 2020-12-25 ezoihy

高斯消元其实就是一个极其靠意识的东西
我们都学过加减消元，在二元时这里其实是极其容易的，但是拓展到多元，我们就需要一种通解
这个东西我在数学课上都听过，我们在考试时也经常用，现在我们用计算机来做其实就更加简单了
以前我们往往是对于一个元两个不同的系数的两个式子，我们往往讲这个元的系数变为原来两个数的系数的最小公倍数，举个例子
[left{ egin{array}{rcl} 5a_1-2a_2=b_12a_1+3a_2=b_2\end{array} ight. ]
我们一般会化成
[left{ egin{array}{rcl} 10a_1-4a_2=2b_110a_1+15a_2=5b_2\end{array} ight. ]

当然我们可以用一般形式来表达一下，我们令
[k*x_i*a_i=y_i*a_i]
显然
[k=frac{y_i}{x_i}]
这样我们将一个式子全部都乘上(k)，我们就可以将(a_i)的系数化成一样的，然后再相减我们就可以消元了
这就是代码，注意在我们判断它的系数极小时将其判定为无解

    for(int i=1;i<=n;i++){
        if(fabs(a[i][i])<1e-7){
            puts("No Solution");
            exit(0);
        }
        const double div=a[i][i];
        for(int j=i;j<=n+1;j++)a[i][j]/=div;
        for(int j=i+1;j<=n;j++){
            const double x=a[j][i];
            for(int k=i;k<=n+1;k++)
                a[j][k]-=a[i][k]*x;
        }
    }

现在我们消完元，我们就可以往回代了，我们已经解出的未知数的值。
其实这里又是一波意识流
代码

    ans[n]=a[n][n+1];
    for(int i=n-1;i;i--){
        ans[i]=a[i][n+1];
        for(int j=i+1;j<=n;j++){
            ans[i]-=a[i][j]*ans[j];
        }
    }

这里我们就完成了消元

luoguP3389 【模板】高斯消元法
完整代码

#include<iostream>
#include<cstdio>
#include<cstring>
#include<cmath>
#include<algorithm>

using namespace std;

int n;
long double ans[101],a[101][101];

void gauss(){
    for(int i=1;i<=n;i++){
        if(fabs(a[i][i])<1e-7){
            puts("No Solution");
            exit(0);
        }
        const double div=a[i][i];
        for(int j=i;j<=n+1;j++)a[i][j]/=div;
        for(int j=i+1;j<=n;j++){
            const double x=a[j][i];
            for(int k=i;k<=n+1;k++)
                a[j][k]-=a[i][k]*x;
        }
    }
    ans[n]=a[n][n+1];
    for(int i=n-1;i;i--){
        ans[i]=a[i][n+1];
        for(int j=i+1;j<=n;j++){
            ans[i]-=a[i][j]*ans[j];
        }
    }
}

int main(){
    scanf("%d",&n);
    for(int i=1;i<=n;i++){
        for(int j=1;j<=n+1;j++){
            double x;
            scanf("%lf",&x);a[i][j]=x;
        }
    }
    gauss();
    for(int i=1;i<=n;i++){
        printf("%.2lf
",(double)ans[i]);
    }
    return 0;
}