高斯消元小结

高斯消元，是一种用于求解线性方程组的算法。

 

对于如图所示的线性方程组，首先我们要将a_i,j与b_i存进一个矩阵，通过初等行变换将其变成容易得出答案的形式。具体实现一般有两种，一种是经典的高斯消元，将其转换成“下三角矩阵”，从末行依次回代得出各未知数的值。另一种是高斯—约旦消元，将所有方程变成kx=b的形式然后直接求解。

高斯—约旦消元（Luogu P3389）

#include<bits/stdc++.h>
using namespace std;
#define il inline
#define ll long long
int read() {
    int s=0,w=1; char ch=getchar();
    while(ch<'0' || ch>'9') { if(ch=='-') w=-1; ch=getchar();}
    while(ch>='0' && ch<='9') s=s*10+ch-'0',ch=getchar();
    return s*w;
} 
const int N=110;
int n;
double a[N][N];
void deal() {
    for(int i=1;i<=n;i++) {
        int maxn=i;
        for(int j=i+1;j<=n;j++) {
            if(fabs(a[j][i])>fabs(a[maxn][i])) maxn=j;
        } //找第i列中最大的系数 
        if(!a[maxn][i]) {
                printf("No Solution\n");
                exit(0);
            } //最大系数为0，无解
        for(int j=1;j<=n+1;j++) swap(a[i][j],a[maxn][j]);//交换第i行与第maxn行
        for(int j=1;j<=n;j++) { 
            if(i==j) continue;
            double temp=a[j][i]/a[i][i];
            for(int k=i+1;k<=n+1;k++) a[j][k]-=a[i][k]*temp;
        } 
    }
}
int main() {
    n=read();
    for(int i=1;i<=n;i++) 
        for(int j=1;j<=n+1;j++) a[i][j]=read();
    deal();
    for(int i=1;i<=n;i++) printf("%.2lf\n",a[i][n+1]/a[i][i]);
    return 0;
}

高斯—约旦消元因为不用回代，所以代码更简洁。但是，这种方法很难处理方程存在自由元（及存在无穷多组解）的情况，所以此时我们会采用更经典的高斯消元。下面的代码处理了无解与存在无穷多解的情况。

高斯消元（Luogu P2455）

#include<bits/stdc++.h>
using namespace std;
#define il inline
#define ll long long
int read() {
    int s=0,w=1; char ch=getchar();
    while(ch<'0' || ch>'9') { if(ch=='-') w=-1; ch=getchar();}
    while(ch>='0' && ch<='9') s=s*10+ch-'0',ch=getchar();
    return s*w;
} 
const int N=310;
const double eps=1e-16;
int n;
double a[N][N];
void deal() {
    for(int i=1;i<=n;i++) {
        int maxn=i;
        for(int j=i+1;j<=n;j++) if(fabs(a[maxn][i])<fabs(a[j][i])) maxn=j;
        if(fabs(a[maxn][i]) < eps) continue; 
        for(int j=1;j<=n+1;j++) swap(a[maxn][j],a[i][j]);
        for(int j=i+1;j<=n;j++) {
            if(i==j) continue;
            double tmp=a[j][i]/a[i][i]; 
            for(int k=i;k<=n+1;k++) {
                a[j][k]-=a[i][k]*tmp;
            }
        }
    }
    for(int i=n;i>=1;i--) {
        if(fabs(a[i][i]) < eps) continue;
        a[i][n+1]/=a[i][i]; 
        for(int j=i-1;j>=1;j--) {
            a[j][n+1]-=a[j][i]*a[i][n+1];
            a[j][i]=0;
        }
    }
}
int main() {
    n=read();
    for(int i=1;i<=n;i++)
        for(int j=1;j<=n+1;j++) a[i][j]=(double)read();
    deal();
    int flag1=0,flag2=0;
    for(int i=1;i<=n;i++) {
        bool flag=true;
        for(int j=1;j<=n;j++) {
            if(fabs(a[i][j]) > eps) {
                flag=false; break;
            }
        }
        if(flag && fabs(a[i][n+1]) > eps) flag1=1;
        if(fabs(a[i][i]) < eps && fabs(a[i][n+1]) < eps) flag2=1;
    }
    if(flag1) {
        printf("-1"); return 0;
    }
    if(flag2) {
        printf("0"); return 0;
    }
    for(int i=1;i<=n;i++) printf("x%d=%.2lf\n",i,a[i][n+1]);
    return 0;
}