高斯消元(求阶梯矩阵,矩阵的秩,解线性方程组)

今天正在学概率dp,在大佬的博客上看到了高斯消元(由于本蒟蒻实在是什么都不懂) 好奇去百度了一下,发现就是手算方程组,然后去洛谷找了个模板题自己手写了一下

题目链接:

https://www.luogu.org/problemnew/show/P3389

很明显高斯消元就是个O(n^3)的暴力算法,这里就不做过多解释了,代码里有少量注释,有线代基础的应该都能看懂

 

下面贴一下自己手写的消元板子(有局限性,因为题里要求没有唯一解就输出No...所以有特别需要的比如求基础解系的可以自行修改)

#include<bits/stdc++.h>

using namespace std;

typedef long long LL;
typedef double D;
D a[105][105];
bool ok = true;
D ans[105];
int n;
template <typename T>
inline void read(T &s){
    T t=1; char k=getchar(); s=0;
    for (;k<'0'||k>'9';k=getchar()) if (k=='-') t=-1;
    for (;k>='0'&&k<='9';k=getchar ()) s=(s<<1)+(s<<3)+(k^48);
    s*=t;
}
inline bool exchange(int j){
    for(int i = j + 1; i <= n; i++){
        if(a[i][j]) swap(a[i], a[j]);
        return true;
    }
    ok = false;//主元个数小于n, 秩小于n 不存在唯一解
    return false;
}
inline void input(){
    for(int i = 1; i <= n;++i){
        for(int j = 1; j <= n + 1; ++j) {
            int tem;
            read (tem);
            a[i][j] = tem;
        }
    }
    D st = a[1][1];
    for(int j = 1; j <= n + 1; ++j){
        a[1][j] /= st;
    }
}
inline void print(){
    for(int i = 1; i <= n; ++i){
        printf("%.2lf\n", ans[i]);
    }
}
void Gauss() {//高斯消元 类比方程组手算
    for (int j = 1; j <= n; ++j) {
        if (a[j][j] || exchange (j)) {//若待消元为0,进行初等行变换
            //要消的是第j个元
            D base = a[j][j];//消元用的基础数
            for (int i = j + 1; i <= n; ++i) {//正在消第i行
                if (a[i][j]) {
                    D multip = a[i][j] / base;//数乘的倍数
                    for (int k = j; k <= n + 1; ++k) {//进行线性运算消元
                        a[i][k] -= a[j][k] * multip;
                    }
                }
            }
        }
        if(!ok) {
            cout << "No Solution"<< endl;
            return;
        }
    }
    //消元完毕 得到一个上三角矩阵
    //进行结果讨论和计算
    for(int i = n; i > 0; --i) {
        for(int j = n; j > i; --j){
            a[i][n + 1] -= a[i][j] * ans[j];
        }
        ans[i] = a[i][n + 1] / a[i][i];
    }
    print();
}
int main(){
    read(n);
    input();
    Gauss();
    return 0;
}

(因为是刚入坑的小萌新,只写过半年代码,所以代码风格或者细节问题求轻喷,代码有问题的话也换迎指正)