【数论学习笔记】高斯消元

模板luogu_3389：https://www.luogu.org/problem/P3389

高斯消元到底是干啥的？？

其实就是解一次方程的。和人一般解方程是类似的，就是让写出来给计算机看就比较麻烦。

还用行列式和矩阵。

先放代码。（照着大佬代码打的……）

#include<bits/stdc++.h>
using namespace std;
double m[111][111];
double ans[111];
double eps=1e-8;
int main(){
    int n;
    cin>>n;
    for(int i=1;i<=n;i++)
        for(int j=1;j<=n+1;j++)
            scanf("%lf",&m[i][j]);
    for(int i=1;i<=n;i++){
        int maxx=i;//找出第一列最大的系数。
        for(int j=i+1;j<=n;j++)
            if(fabs(m[maxx][i])<fabs(m[j][i]))
                maxx=j;
        if(fabs(m[maxx][i])<eps){
            printf("No Solution");//无解。
            return 0;
        }
        if(i!=maxx) 
            swap(m[i],m[maxx]); //交换第一行和系数最大行。
        double div=m[i][i];//化简这一横行系数，让除0外第一个为1。
        for(int j=i;j<=n+1;j++)
            m[i][j]/=div;
        for(int j=i+1;j<=n;j++){
            div=m[j][i];
            for(int k=i;k<=n+1;k++)
                m[j][k]-=m[i][k]*div;
        }
    }
    ans[n]=m[n][n+1];
    for(int i=n-1;i>=1;i--){
        ans[i]=m[i][n+1];
        for(int j=i+1;j<=n;j++)
            ans[i]-=(m[i][j]*ans[j]);
    }
    for(int i=1;i<=n;i++)
        printf("%.2lf
",ans[i]);
    // while(1);
    return 0;
}

具体解析就不打了，在蓝书P158页。还有洛谷某大佬的题解：https://pks-loving.blog.luogu.org/p3389-mu-ban-gao-si-xiao-yuan-fa

但最后的ans数组比较重要，题解大佬说这是倒带操作。

实际上就是把之前没有减的常数在最后减了。

因为要对角线，所以ans[i]一定要减了后面行的后面列的系数倍。

emmmm，不太清楚，但是也理解一下。

UPDATE 19.8.20

感性理解了很长时间，突然记起来不清楚的时候还写了题解。

在今天复习心头一热，码一个模板题目。

最后那个倒带其实不是什么高深操作。

就是在消元形成上三角的时候，除去之前算出的答案。

就是这一行的未知数的答案。

感性理解，感性理解……

19.11.4

开始用高斯约旦消元了

还挺好写的

代码是LuoguP2455 [SDOI2006]线性方程组

代码如下:

#include<bits/stdc++.h>
using namespace std;
const int maxn=100;
const double eps=1e-8;
int n,id[maxn];
double a[maxn][maxn];
int main()
{
    scanf("%d",&n);
    for(int i=1;i<=n;i++) for(int j=1;j<=n+1;j++)
        scanf("%lf",&a[i][j]);
    int pos=1;
    for(int i=1;i<=n;i++){
        int maxx=pos;
        for(int j=pos+1;j<=n;j++) if(fabs(a[j][i])>fabs(a[maxx][i])) maxx=j;
        if(!a[maxx][i]) continue;
        swap(a[pos],a[maxx]);
        for(int j=1;j<=n;j++){
            if(j==pos) continue;
            double t=1.0*a[j][i]/a[pos][i];
            for(int k=i;k<=n+1;k++) a[j][k]-=t*a[pos][k];
        }
        id[i]=pos;pos++;
    }
    for(int i=1;i<=n;i++) if(!id[i]) id[i]=(pos==n) ? 0 : pos++;
    for(int i=1;i<=n;i++) if(!a[id[i]][i] && a[id[i]][n+1]){puts("-1");return 0;}
    for(int i=1;i<=n;i++) if(!a[id[i]][i] && !a[id[i]][n+1]){puts("0");return 0;}
    for(int i=1;i<=n;i++)
        printf("x%d=%.2lf
",i,fabs(a[id[i]][n+1]/a[id[i]][i])<eps ? 0.0 : a[id[i]][n+1]/a[id[i]][i]);
    return 0;
}