2018-01-24

高斯消元小结

介绍

高斯消元是一种解线性方程组的 $O(N ^ 3)$ 的方法。其实也就是线性代数中说的将一个矩阵化为上三角矩阵的的过程，最后再求解每个未知量的值。

类比下，高斯消元的在解方程组的作用，和最短路在查分约束中的作用、网络流在网络流的图论模型中的作用差不多。都只是一个计算的工具，题目考察的核心点并不在高斯消元上而是在如何列出对应的线性方程组。

另外要说的高斯消元大致有三种写法，第一种是简单的，矩阵的值只有 0 和 1，这里高斯消元中只用异或运算就可以了，其实就是求解线性基的过程；第二种是解一些模线性方程组，答案需要整数解，这里消元的时候需要取下最小公倍数；最后一种是没有特殊要求，直接上 Double 类型的，解决概率期望的问题会用到。

消元过程

01矩阵

int solve() {
    //rst 代表确定的未知量个数，也就是矩阵的秩、最大无关组
    //r 代表当前进行到了第几行，之前的行都化为了阶梯型
    int rst = 0, r = 0;
    memset(lib, 0, sizeof(lib));
    //代表当前进行到了第几列
    for (int i = 0; i < n; i++) {
        //枚举这一列的所有行，将一个 i 列不为 0 的行换到第 r 行
        for (int j = r; j < n; j++) {
            if (!data[j][i]) continue;
            for (int k = i; k <= n; k++) {
                swap(data[j][k], data[r][k]);
            }
            break;
        }
        //如果当前行的第 i 列全都为 0 那么就跳过这一列
        if (!data[r][i]) continue;
        rst++;
        //用第 r 行去消去其他行，如果其他行第 i 列是 1 的话
        for (int j = 0; j < n; j++) {
            if (j == r || !data[j][i]) continue;
            for (int k = i; k <= n; k++) data[j][k] ^= data[r][k];
        }
        lib[i] = true, r++;
    }
}

模线性方程组

前面的部分都一样，只是进行消元的地方修改一下，求个最小公倍数。

for (int j = 0; j < n; j++) {
    if (j == r || !data[j][i]) continue;
    int t = lcm(data[r][i], data[j][i]);
    int t1 = t / data[r][i];
    int t2 = t / data[j][i];
    for (int k = 0; k <= n; k++) {
        data[j][k] = ((data[j][k] * t2 - data[r][k] * t1) % MOD + MOD) % MOD;
    }
}

浮点数矩阵

消元的地方直接求个浮点数的倍数就可以。

for (int j = 0; j < n; j++) {
    if (j == r || eq(data[j][i], 0)) continue;
    double temp = data[j][i] / data[r][i];
    for (int k = i; k <= n; k++) {
        data[j][k] -= data[r][k] * temp;
    }
}

处理答案

高斯消元算法执行结束后，会得到如下信息：自由未知量的个数，那些列是自由未知量，方程的一个解。

01矩阵处理答案的过程是这样的，这里默认让所有自由未知量为 0 。

for (int i = 0; i < n; i++) {
    //如果第 n 列不为 0，其他列全为 0，说明增广矩阵的秩大于原矩阵，方程无解
    if (data[i][n]) {
        bool flag = false;
        for (int j = 0; j < n; j++) {
            if (data[i][j]) flag = true;
        }
        if (!flag) return -1;
    }
    //如果有自由未知量直接跳过
    if (!l[i]) continue;
    //否则这个变量的值就等于第 n 列的值，因为当前这一列除了自由未知量的列，其他列全为 0
    for (int j = 0; j < n; j++) {
        if (data[j][i]) ans[i] = data[j][n];
    }
}
return n - rst;

01矩阵有时候需要的一种做法是枚举所有的自由未知量，求取 1 的变量最少的解。

一般是直接用二进制去枚举。

int anss = INF;
for (int k = 0; k < 16; k++) {
    int cnt1 = 0, cnt2 = 0;
    int num = (k & 1) + (k >> 1 & 1) + (k >> 2 & 1) + (k >> 3 & 1);
    for (int i = 0; i < 12; i++) {
        int t = 0;
        for (int j = 0; j < 4; j++) if (data[i][12 + j]) t ^= ((k >> j) & 1);
        if (t != data[i][n]) cnt1++;
        if (t != data[i][n + 1]) cnt2++;
    }
    anss = min(min(cnt1, cnt2) + num, anss);
}

Fish's Blog