[线性代数]高斯消元(Gaussian Elimination)

在算法竞赛中，只要遇到“把一个向量v变成另一个向量v’，并且v’的每一个分量都是v各个分量的线性组合”的情况，就可以考虑用矩阵乘法来描述这个关系。
高斯消元可以用来求行列式、线性方程组、异或方程组、逆矩阵、伴随矩阵、矩阵的秩、自由变元个数，是个很强的数学工具。

高斯消元求解线性方程组(混合高斯约当消元法)(极速版)：
https://www.luogu.org/problemnew/show/P3389

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int maxn = 4e2 + 5, mod = 1e9 + 7;
const double eps = 1e-10;
double a[maxn][2 * maxn];

bool gauss(int n, int m)
{
    for (int i = 0; i < n; ++i){

        int r = i;
        for (int j = i + 1; j < n; ++j){
            if (fabs(a[j][i]) > fabs(a[r][i])) r = j;
        }
        if (r != i){
            for (int j = i; j < m; ++j) swap(a[r][j], a[i][j]);
        }
       
        if (a[i][i] == 0) return false;
        for (int k = i + 1; k < n; ++k){
            double f2 = a[k][i] / a[i][i];
            for (int j = i; j < m; ++j){
                a[k][j] -= f2 * a[i][j];
            }
        }
        /*
        //更精确的写法
        for (int j = m - 1; j >= i; --j){
            for (int k = i + 1; k < n; ++k){
                a[k][j] -= a[k][i] / a[i][i] * a[i][j];
            }
        }*/
    }

    return true;
}
//这就是把系数矩阵直接变成单位矩阵
bool GJ(int n, int m)
{
    if (!gauss(n, m)) return false;
    for (int i = 0; i < n; ++i){
        if (a[i][i] == 0) return false;
        for (int k = 0; k < i; ++k){
            double f2 = a[k][i] / a[i][i];
            for (int j = i; j < m; ++j){
                a[k][j] -= f2 * a[i][j];
            }
        }
    }
    for (int i = 0; i < n; ++i){
        for (int j = n; j < m; ++j) a[i][j] = a[i][j] / a[i][i];
    }
    return true;
}
int main()
{
    ios::sync_with_stdio(0); cin.tie(0);
    int n;
    cin >> n;
    for (int i = 0; i < n; ++i){
        for (int j = 0; j < n + 1; ++j){
            cin >> a[i][j];
        }
    }
    if (!GJ(n, n + 1)) cout << "No Solution" << '\n';
    else{
        for (int i = 0; i < n; ++i) cout << fixed << setprecision(2) << a[i][n] << '\n';
    }
    return 0;
}

高斯消元求逆矩阵对1e9+7取模(极速版)：
https://www.luogu.org/problemnew/show/P4783

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int maxn = 4e2 + 5, mod = 1e9 + 7;
const double eps = 1e-10;
int a[maxn][2 * maxn];

ll quickpow(ll x, int k)
{
    ll res = 1;
    while (k){
        if (k & 1) res = (res * x) % mod;
        x = (x * x) % mod;
        k >>= 1;
    }
    return res;
}
bool gauss(int n, int m)
{
    for (int i = 0; i < n; ++i){

        int r = i;
        for (int j = i + 1; j < n; ++j){
            if (abs(a[j][i]) > abs(a[r][i])) r = j;
        }
        if (r != i){
            for (int j = i; j < m; ++j) swap(a[r][j], a[i][j]);
        }

        if (a[i][i] == 0) return false;
        ll f1 = quickpow(a[i][i], mod - 2);
        for (int k = i + 1; k < n; ++k){
            ll f2 = a[k][i] * f1 % mod;
            for (int j = i; j < m; ++j){
                a[k][j] -= f2 * a[i][j] % mod;
                a[k][j] = (a[k][j] + mod) % mod;
            }
        }

    }

    return true;
}
//这就是把系数矩阵直接变成单位矩阵
bool GJ(int n, int m)
{
    if (!gauss(n, m)) return false;
    for (int i = 0; i < n; ++i){
        if (a[i][i] == 0) return false;
        ll f1 = quickpow(a[i][i], mod - 2);
        for (int k = 0; k < i; ++k){
            ll f2 = a[k][i] * f1 % mod;
            for (int j = i; j < m; ++j){
                a[k][j] -= f2 * a[i][j] % mod;
                a[k][j] = (a[k][j] + mod) % mod;
            }
        }
    }
    for (int i = 0; i < n; ++i){
        ll f1 = quickpow(a[i][i], mod - 2);
        for (int j = n; j < m; ++j) a[i][j] = a[i][j] * f1 % mod;
    }
    return true;
}
int main()
{
    ios::sync_with_stdio(0); cin.tie(0);
    int n;
    cin >> n;
    for (int i = 0; i < n; ++i){
        for (int j = 0; j < n; ++j){
            cin >> a[i][j];
        }
    }
    for (int i = 0; i < n; ++i) a[i][i + n] = 1;
    if (!GJ(n, 2 * n)) cout << "No Solution" << '\n';
    else{
        for (int i = 0; i < n; ++i){
            for (int j = n; j < 2 * n; ++j){
                cout << a[i][j] << ' ';
            }
            cout << '\n';
        }
    }
    return 0;
}

高斯消元解异或方程组求自由变元个数：
http://acm.hdu.edu.cn/showproblem.php?pid=5833
自由变元个数 = 自变量个数 – 矩阵的秩

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int maxn = 2e3 + 5, mod = 1e9 + 7;

ll ta[305];
int a[maxn][305], seive[maxn], sp;
bool prime[100005];

void print(int n, int m)
{
    cout << endl;
    for (int i = 0; i < n; ++i){
        for (int j = 0; j < m; ++j){
            cout << a[i][j] << ' ';
        }
        cout << endl;
    }
}
void init()
{
    for (int i = 2; i < 100005; ++i){
        if (!prime[i]){
            seive[sp++] = i;
            if (sp > 2000) break;
            for (int j = 2; i * j < 100005; ++j){
                prime[i * j] = 1;
            }
        }
    }
}
int dis(int n)
{
    int p = 0;
    for (int i = 0; i < n; ++i){
        for (int j = 0; j < sp; ++j){
            if (seive[j] * seive[j] > ta[i]) break;
            int cnt = 0;
            while (ta[i] % seive[j] == 0){
                cnt ^= 1;
                ta[i] /= seive[j];
            }
            if (cnt) p = max(p, j);
            a[j][i] = cnt;
        }
        if (ta[i] != 1){
            int t = lower_bound(seive, seive + sp, ta[i]) - seive;
            p = max(p, t);
            a[t][i] = 1;
        }
    }
    return p;
}
int quickpow(ll x, int k)
{
    ll res = 1;
    while (k){
        if (k & 1) res = (res * x) % mod;
        x = (x * x) % mod;
        k >>= 1;
    }
    return res;
}
int gauss(int n, int m)
{
    int p = 0;//当前应处理的行数
    //当前处理第i列
    for (int i = 0; i < m; ++i){
    //从第p行开始找
        int r = p;
        while (r < n){
            if (a[r][i] == 1) break;
            else ++r;
        }
        //找到为1的，直接提上来
        if (a[r][i] == 0) continue;
        if (r != p) for (int j = i; j < m; ++j) swap(a[r][j], a[p][j]);
        //把下面为1的整条异或掉
        for (int k = p + 1; k < n; ++k){
            if (a[k][i] == 1){
                for (int j = i; j < m; ++j){
                    a[k][j] ^= a[p][j];
                }
            }
        }
        //p行已变换完毕
        ++p;
    }
    return p;
}

int main()
{
    ios::sync_with_stdio(0); cin.tie(0);
    init();
    int t;
    cin >> t;
    for (int cas = 1; cas <= t; ++cas){
        int n;
        cin >> n;
        for (int i = 0; i < n; ++i) cin >> ta[i];
        int num = dis(n);
        int p = n - gauss(num + 1, n);
        cout << "Case #" << cas << ":\n" << quickpow(2, p) - 1 << '\n';
        memset(a, 0, sizeof(a));
    }
    return 0;
}

生成树计数：
对于一个无向图 G ，它的生成树个数等于其基尔霍夫Kirchhoff矩阵任何一个N-1阶主子式的行列式的绝对值。
基尔霍夫Kirchhoff矩阵 K =度数矩阵 D – 邻接矩阵 A
求行列式高斯消元即可