高斯消元
upd：2024/7/14 几乎是重构本文章。
upd：2024/7/15 修改了本文章的部分描述。
前言：参考自《高等代数（上册）》。本文适用于中国宝宝体质食用。
一、算法介绍
定义1： 一个 \(n \times m\) 的矩阵 \(A\) 如果满足以下条件就被称为阶梯形矩阵。
- 零行（元素全为 \(0\) 的行）在最下方。
- 主元（非零行的第一个不为 \(0\) 的元素）的列随行的增大而增大。
定义2： 将以下三种变换称为矩阵的初等行变换。
- 1° 型：某行加上另外一行的某个倍数
- 2° 型：交换任意两行
- 3° 型：用一个非零数乘上某一行
命题1： 将初等行变换把线性方程组变成同解的方程组。以 1° 型为例，设原方程组为 \(A\)，变换一次后为 \(B\)，容易证明 \(A\) 的每一个解都是 \(B\) 的解，将 \(B\) 变换回去又可以证明 \(B\) 的每一个解都是 \(A\) 的解，因此 \(A\) 与 \(B\) 同解。
命题2： 任意一个矩阵都可以通过初等行变换化成阶梯形矩阵。这个可以通过数学归纳法证明。
当我们把方程组（设有 \(m\) 个方程，\(n\) 个未知数）化为阶梯形矩阵后，就不难得出解了。如果存在 \(0=d \ (d \ne 0)\) 的行，原方程组无解。否则设 \(r\) 为阶梯形矩阵中不为 \(0=0\) 型的行的个数。分为以下两种情况：
- \(r=n\) ：此时原方程组有唯一解，且前 \(n\) 行一定为上三角的形式，从下往上依次带入计算即可。
- \(r<n\) ：我们称不为 \(0=0\) 型的行的主元为系数的未知量为主变量，其余未知量称为自由未知量。可以发现对自由未知量任意赋一组值后，所有主变量的值也就唯一确定了。因此这样的方程组具有无穷个解。
除了通过将矩阵变换为阶梯形矩阵来求解外，还可以变换为简化阶梯形矩阵。
定义3：如果一个阶梯形矩阵满足其所有主元全为 \(1\)，且每个主元所在的列的其余元素均为 \(0\)，那么称它为简化阶梯形矩阵。
这样做的优点是当 \(r=n\) 时不需要回带。以上便是高斯消元的两种方法。该算法的时间复杂度为 \(O(n^2m)\) 的。
以下就是 P2455 [SDOI2006] 线性方程组 的代码（注释里附有实现细节）：
#include<bits/stdc++.h>
#define int long long
#define N 110
#define lb long double
using namespace std;
int n;
lb a[N][N], b[N]; //a[i][j]表示第i个方程,第j个系数
lb x[N];
lb Abs(lb x) {
return x > 0 ? x : -x;
}
void My_work() { //第二种消元，也被称为约旦消元QAQ
int r = 1; //r表示0=0和0=d(d!=0)的行的个数（具体实现时考虑到有无解的情况与文章定义略有不同）。
for(int i = 1; i <= n; i++) { //每个未知数
lb Mx = 0;
int id = 0;
for(int j = r; j <= n; j++) { //第j个方程
if(Abs(a[j][i]) > Mx) { //选一个绝对值最大的（据说这样做能降低精度）
Mx = a[j][i];
id = j;
}
}
if(id == 0) continue;
for(int j = i; j <= n + 1; j++) swap(a[r][j], a[id][j]); //把找到的行放在第r行
for(int j = 1; j <= n; j++) { //将第j个方程的第i个方程的系数消成0
if(a[j][i] == 0 || j == r) continue;
lb Num = 1.0 * a[j][i] / a[r][i];
for(int k = 1; k <= n + 1; k++) a[j][k] -= a[r][k] * Num;
}
/*for(int i = 1; i <= n; i++) {
for(int j = 1; j <= n + 1; j++) printf("%.2Lf ", a[i][j]);
printf("\n");
}
printf("\n");*/
r++;
}r--;
//printf("%lld\n",r);
if(r < n) {
for(int i = r + 1; i <= n; i++) {
if(Abs(a[i][n + 1]) > 1e-9) { //无解
cout << -1 << endl;
return;
}
}
cout << 0 << endl; //无数解
return;
}
for(int i = 1; i <= n; i++)
printf("x%lld=%.2Lf \n", i, a[i][n + 1] / a[i][i]);
}
signed main() {
//ios::sync_with_stdio(0);
//cin.tie(0); cout.tie(0);
scanf("%lld", &n);
for(int i = 1; i <= n; i++) {
for(int j = 1; j <= n + 1; j++) scanf("%Lf", &a[i][j]);
}
My_work();
return 0;
}
二、例题讲解
1、P3164 [CQOI2014] 和谐矩阵
题意：
我们称一个由 \(0\) 和 \(1\) 组成的矩阵是和谐的，当且仅当每个元素都有偶数个相邻的 \(1\)。一个元素相邻的元素包括它本身，及他上下左右的 \(4\) 个元素（如果存在）。给定矩阵的行数和列数，请计算并输出一个和谐的矩阵。注意：所有元素为 \(0\) 的矩阵是不允许的。
\(1\le n,m\le 40\)。
分析：
构造异或线性方程组，对高斯消元稍加修改就能支持了。利用 bitset 优化，时间复杂度 \(O(\frac{n^2m^2}{w})\)，但暴力消元也能通过本题。
2、P2447 [SDOI2010] 外星千足虫
题意：
现在你面前摆有 \(1\ldots N\) 编号的 \(N\) 只千足虫，你的任务是鉴定每只虫子所属的星球，但不允许亲自去数它们的足。
Charles 每次会在这 \(N\) 只千足虫中选定若干只放入“昆虫点足机”（the Insect Feet Counter, IFC）中，“点足机”会自动统计出其内所有昆虫足数之和。Charles 会将这个和数 \(\bmod\) \(2\) 的结果反馈给你，同时告诉你一开始放入机器中的是哪几只虫子。
他的这种统计操作总共进行 \(M\) 次，而你应当尽早得出鉴定结果。
假如在第 \(K\) 次统计结束后，现有数据就足以确定每只虫子的身份，你就还应将这个 \(K\) 反馈给 Charles，此时若 \(K<M\)，则表明那后 \(M-K\) 次统计并非必须的。
如果根据所有 \(M\) 次统计数据还是无法确定每只虫子身份，你也要跟 Charles 讲明：就目前数据会存在多个解。
对于 \(100\%\) 的数据，满足 \(1\leq N\leq 10^3\)，\(1\leq M\leq 2\times 10^3\)。
分析：
\(\mod2\) 这个条件同样可以转异或方程组。保证使用的方程最少，只需要选主元时贪心取编号最小的方程即可。同样也能用 bitset 优化。
