高斯消元 学习笔记

发表于 | 分类于 | | + s

模板及讲解

高斯消元是什么

求解多元一次方程组的一个消元方法，主要借助矩阵来完成。

例题

Luogu 2455

我们发现给定的方程组解只和系数有关，所以将系数和常数提取出来构造矩阵：
例
$$
\begin{cases}
2x+y+z=1\\
6x+2y+z=-1\\
-2x+2y+z=7
\end{cases}
$$
则构造矩阵：
$$
\left[
\begin{array}{ccc|c}
2 & 1 & 1&1 \\
6 & 2 & 1&-1\\
-2&2&1&7
\end{array}
\right]
$$

高斯消元的操作有(初等行变换)：
1、交换两行。
2、将一行用非零数乘。
3、将一行的若干倍加到另一行上。

也就是差不多是小学学的加减消元带入消元的操作。

高斯消元目标：
将原矩阵转化为一个阶梯形矩阵，即
$$
\left[
\begin{array}{ccc|c}
2 & 1 & 1&1 \\
0 & -1& -2&-4\\
0&0&-4&-4
\end{array}
\right]
$$
将后面的行带入上面的行，可得一个简化阶梯形矩阵
$$
\left[
\begin{array}{ccc|c}
1 & 0 & 0&-1 \\
0 & 1& 0&2\\
0&0&1&1
\end{array}
\right]
$$
直接可以读出每个元素的解。

高斯消元一般步骤：
逐步处理每一个$x_i$，将$i$列除了$i$行的其他系数都变为0，即我们可以先保证$i$行$i$列不为0，然后以这一行为基准将其他行$i$位置的系数都变为0。用初等行变换的交换可以满足第一个条件，用初等行变换将一行的若干倍加到另一行上可以做到第二个条件。做完之后只有对角线和常数有值，每一行的常数除以对角线的值即为第$i$个元素的解。

无解情况：系数都为0，常数不为0
无穷解情况：系数都为0，常数也为0

两个概念：
1、主元：确定唯一的未知量
2、自由元：不确定值未知量，对应情况为该行全为0，有无穷解

注意在判这两个之前要先判有没有不是对角线上的系数有值，如果有则这一行不进行检查。

#include<cstdio> 
#include<cstring>
#include<algorithm>
#include<iostream>
#include<vector>
#include<queue>
#include<cmath>
#define ms(i, j) memset(i, j, sizeof i)
#define LL long long
#define db double
#define fir first
#define sec second
#define mp make_pair
using namespace std;

namespace flyinthesky {

    const int MAXN = 105;
    const db eps = 1e-12; 

    int n;
    db c[MAXN][MAXN], b[MAXN];
    // c, b 为消元矩阵 

    void clean() {
        ms(b, 0);
    }
    int solve() {
        clean();
        scanf("%d", &n);
        for (int i = 1; i <= n; ++i) {
            for (int j = 1; j <= n; ++j) scanf("%lf", &c[i][j]);
            scanf("%lf", &b[i]);
        }
        for (int i = 1; i <= n; ++i) { // 处理 xi[i] 的系数 
            for (int j = i; j <= n; ++j) { // 从后面没处理的行找 xi[i] 系数不为 0 的行交换 
                if (fabs(c[j][i]) > eps) {
                    swap(b[i], b[j]);
                    for (int k = 1; k <= n; ++k) swap(c[i][k], c[j][k]);
                }
            }
            if (fabs(c[i][i]) < eps) continue ; // 注意 
            for (int j = 1; j <= n; ++j) { // 消去其他方程 xi[i] 的系数 
                if (i == j) continue ;
                db rat = c[j][i] / c[i][i];
                for (int k = 1; k <= n; ++k) c[j][k] -= rat * c[i][k];
                b[j] -= rat * b[i];
            }
        }
        int bj = 1;
        for (int i = 1; i <= n; ++i) {
            int fl = 0;
            for(int j = 1; j <= n; ++j){
                if (i != j && fabs(c[i][j]) > eps) {fl = 1; break ;}
            }
            if (fl) continue ; // 必须判 
            if (fabs(c[i][i]) < eps && fabs(b[i]) > eps) bj = -1;
            if (fabs(c[i][i]) < eps && fabs(b[i]) < eps && bj != -1) bj = 0; 
        }
        if (bj != 1) return printf("%d\n", bj), 0; 
        for (int i = 1; i <= n; ++i) {
            if (fabs(b[i] / c[i][i]) < eps) printf("x%d=0\n", i); else printf("x%d=%.2lf\n", i, b[i] / c[i][i]);
        }
        return 0;
    }
}
int main() {
    flyinthesky::solve();
    return 0;
}

异或方程组

高斯消元可以用来解线性的异或方程组，方法与前者类似。

例题：poj 1830

设$x_i$为$i$开关是否(0/1)操作，$a_{i,j}$为$j$开关操作后是否(1/0)会影响$i$开关

则可以列出线性异或方程组

$$
\begin{cases}
a_{1,1}x_1 xor a_{1,2}x_2 xor … xor a_{1,n}x_n = st_1 xor ed_1 \\
a_{2,1}x_1 xor a_{2,2}x_2 xor … xor a_{2,n}x_n = st_2 xor ed_2 \\
\vdots \\
a_{n,1}x_1 xor a_{n,2}x_2 xor … xor a_{n,n}x_n = st_n xor ed_n
\end{cases}
$$

然后对于异或实际上我们也可以进行消元得到简化阶梯形矩阵，并且不需要乘除，因为只有 0 和 1。
可以用 bitset 加速运算且代码更简洁高效。

无解即为形如$0 xor 0 xor … xor 0 =1$的，显然不成立。

#include<cstdio> 
#include<cstring>
#include<algorithm>
#include<iostream>
#include<vector>
#include<queue>
#include<bitset>
#define ms(i, j) memset(i, j, sizeof i)
#define LL long long
#define db double
#define fir first
#define sec second
#define mp make_pair
using namespace std;

namespace flyinthesky { 

    const int MAXN = 35;

    int n, a[MAXN];
    bitset<MAXN> bs[MAXN];

    int ksm(int a, int b) {
        int bs = a, ans = 1;
        while (b) {
            if (b & 1) ans *= bs;
            bs *= bs;
            b >>= 1;
        }
        return ans;
    }

    void clean() {
    }
    int solve() {
        clean();
        cin >> n;
        for (int i = 1; i <= n; ++i) scanf("%d", &a[i]);
        for (int x, i = 1; i <= n; ++i) scanf("%d", &x), a[i] ^= x;
        for (int i = 1; i <= n; ++i) bs[i].reset(), bs[i][0] = a[i], bs[i][i] = 1;
        int x, y;
        while (scanf("%d%d", &x, &y) == 2 && x && y) bs[y][x] = 1;
        int ans = 1, whw = 0;
        for (int i = 1; i <= n; ++i) { // 每列
             if (bs[i][i] == 0) { // 找一个非零的
                 for (int j = i + 1; j <= n; ++j) if (bs[j][i]) {swap(bs[j], bs[i]); break;}
             }
             if (bs[i][0] == 1 && bs[i].count() == 1) {ans = 0; break ;} // 无解
             if (bs[i][i] == 0) ++whw; // 自由元
             for (int j = i + 1; j <= n; ++j) if (bs[j][i]) bs[j] ^= bs[i];
        }
        if (ans == 0) printf("Oh,it's impossible~!!\n"); else printf("%d\n", ksm(2, whw));
        return 0; 
    }
}
int main() {
    int k; cin >> k;
    while (k--) flyinthesky::solve();
    return 0;
}

同余高斯消元

模板题：CF 1155E

求解带环期望DP

1、CF 24D：列与列之间不满足无后效性 (注意此时高斯消元做成阶梯矩阵即可，不需要简化，简化直接回带)

线性基

定义

一下内容部分来自 线性基学习笔记 | Menci’s Blog

异或和：设 $S$ 为一无符号整数集（下文中除非特殊说明，集合均指无符号整数集），定义其异或和为$S_1 xor S_2 xor \cdots xor S_{|S|}$。

线性相关：集合中存在一个元素$S_j$，使得集合中其他元素可以通过异或和的方式得到$S_j$，则$S$线性相关。若不存在，则$S$线性无关。

线性基：称集合 $B$ 是集合 $S$ 的线性基，当且仅当：

  • $B$ 线性无关
  • $B$ 中元素可以通过异或和方式得到$S$中的所有值

线性基性质：

  • $B$ 是极小的满足线性基性质的集合，它的任何真子集都不可能是线性基
  • $S$ 中的任意元素都可以唯一表示为 $B$ 中若干个元素异或起来的结果。

构造

构造出来的线性基相当于于所有插入值高斯消元后的简化阶梯矩阵。

设维护的数的二进制长度为$L$，线性基数组为$a_i(0 \leq i \leq L)$
为了保证线性无关，$a_i$的第$i$位若为$1$，则其他线性基上第$i$位都为$0$，并且$a_i$高于$i$位上的值都必须为$0$(否则会影响后面线性基的”其他线性基上第$i$位都为$0$”)。

所以构造方法即为(插入数$t$)：
1、将$t$倒序枚举二进制位$i$，然后
2、若$t$的$i$位为0，则continue
否则：

  • 如果$a_i$不为$0$，则用$a_i$消去$t$的$i$位上的二进制$1$，然后continue
  • 如果$a_i$为$0$，那么$t$就能插进$i$位。但是要先用$[1, i)$的线性基消去$t$的相应位置的$1$，以保证”其他线性基上第$i$位都为$0$”，然后再用得到的线性基消去$[i+1, L]$线性基的$i$位上的$1$，因为一路上我们都将高位消去，所以这里不用担心影响其他位的情况。

例题1 (集合异或和最大)：Luogu 3812

给定$n$个整数（数字可能重复），求在这些数中选取任意个，使得他们的异或和最大。

对于构造的线性基，$a_i$的第$i$位若是$1$，则将答案异或$a_i$，这样一定能得到更大的值，因为这个位置的$1$之后不会再被异或消去，而因为二进制的性质，贪心选取更高位数的$1$，所以将所有线性基异或起来即为答案。

#include<cstdio> 
#include<cstring>
#include<algorithm>
#include<iostream>
#include<vector>
#include<set>
#include<cmath>
#define ms(i, j) memset(i, j, sizeof i)
#define LL long long
#define db double
#define fir first
#define sec second
#define mp make_pair
using namespace std;

namespace flyinthesky {

    const LL MAXLEN = 51;

    LL n, t, a[55];

    void clean() {
        ms(a, 0);
    }
    int solve() {

        clean();

        scanf("%lld", &n);

        for (LL o = 1; o <= n; ++o) {
            scanf("%lld", &t);
            for (LL i = MAXLEN; i >= 0; --i) { // 倒序枚举二进制
                if (!((t >> i) & 1)) continue ; // t 当前位不是 1
                if (a[i]) t ^= a[i]; // 用 a[i] 消去 t 的 i 位上的 1
                  else { 
                    for (LL j = 1; j < i; ++j) if ((t >> j) & 1) t ^= a[j]; // 用 a[j] 消去 t 低于 i 位上的 1
                    for (LL j = i + 1; j <= MAXLEN; ++j) if ((a[j] >> i) & 1) a[j] ^= t; // 用 t 消去 a[j] 的 i 位上的 1
                    a[i] = t; // 插入线性基
                    break ;
                }
            }
        }
        LL ans = 0;
        for (LL i = MAXLEN; i >= 0; --i) ans ^= a[i];

        printf("%lld\n", ans);

        return 0;
    }
}
int main() {
    flyinthesky::solve();
    return 0;
}

例题2 (选出集合不出现异或和等于0)：Bzoj 3105

选出集合不出现异或和等于0即选出数线性无关，所以直接最大到最小贪心即可。

#include<cstdio> 
#include<cstring>
#include<algorithm>
#include<iostream>
#include<cmath>
#include<map>
#include<queue>
#include<string>
#define ms(i, j) memset(i, j, sizeof i)
#define LL long long
#define db double
#define fir first
#define sec second
#define mp make_pair
using namespace std;

namespace flyinthesky {

    int k, a[105], b[105];

    int ins(int t) {
        for (int i = 32; ~i; --i) {
            if (!(t & (1 << i))) continue ;
            if (a[i]) t ^= a[i];
            else {
                for (int j = 0; j < i; ++j) if (t & (1 << j)) t ^= a[j];
                for (int j = i + 1; j <= 32; ++j) if (a[j] & (1 << j)) a[j] ^= t;
                a[i] = t;
                return 1;
            }
        }
        return 0;
    }

    void clean() {
    }
    int solve() {

        clean();
        cin >> k;
        for (int i = 1; i <= k; ++i) cin >> b[i];
        sort(b + 1, b + 1 + k, greater<int>());

        LL ans = 0;

        for (int i = 1; i <= k; ++i) if (!ins(b[i])) ans += b[i];

        cout << ans;

        return 0;
    }  
}
int main() {
    flyinthesky::solve();
    return 0;
}

例题3 求集合异或和$k$大

我们考虑将秩编号，那么我们如果求$k$大, 我们就分解$k$的二进制, 然后按照二进制的位置异或秩就是答案，注意0的存在

例题4 求线性基的交 / 并

并很简单，一个插另一个即可

交

例题5 求集合异或和排名

类似例题三的做法，Bzoj 2844

此题还需要去重，我们发现每个数都出现一样的次数$2^{n-|\mathbb{lb}|}$

证明：假如要考虑异或出一个数$x$，基外选出的数的异或和为$s$，那么还需要$x \mathbb{xor} s$，对于它，基内的线性表示的方法是唯一的。

------ 本文结束 ------