






















在算法竞赛中,高斯消元法常用于求解 \(n\) 元一次方程组。
一个包含 \(n\) 个变量和 \(m\) 个方程的线性方程组可以表示为:
\[\begin{cases} a_{11}x_1 + a_{12}x_2 + \dots + a_{1n}x_n = b_1 \\ a_{21}x_1 + a_{22}x_2 + \dots + a_{2n}x_n = b_2 \\ \dots \\ a_{m1}x_1 + a_{m2}x_2 + \dots + a_{mn}x_n = b_m \end{cases} \]
其对应的增广矩阵为:
\[\begin{bmatrix} a_{11} & a_{12} & \dots & a_{1n} & b_1 \\ a_{21} & a_{22} & \dots & a_{2n} & b_2 \\ \vdots & \vdots & \ddots & \vdots & \vdots \\ a_{m1} & a_{m2} & \dots & a_{mn} & b_m \end{bmatrix} \]
利用以下三种变换(初等行变换),可以将增广矩阵转化为行阶梯形矩阵或简化行阶梯形矩阵:
通过这些变换,将复杂的方程组转化为一个上三角形(行阶梯形)或对角形(简化行阶梯形)的方程组。例如,将方程组转化为简化行阶梯形矩阵后:
\[\begin{bmatrix} 1 & 0 & \dots & 0 & c_1 \\ 0 & 1 & \dots & 0 & c_2 \\ \vdots & \vdots & \ddots & \vdots & \vdots \\ 0 & 0 & \dots & 1 & c_n \end{bmatrix} \]
这直接对应了方程组的解:
\[\begin{cases} x_1 = c_1 \\ x_2 = c_2 \\ \dots \\ x_n = c_n \end{cases} \]
因此,高斯消元的过程就是通过一系列保解的变换,将方程组简化为可以直接读出解的形式。
以下是 Gauss-Jordan 消元法的算法流程:
矩阵的秩就等于高斯消元后阶梯形矩阵中非零行的行数。高斯消元的本质是通过行与行的加减消去多余的信息,如果某一行可以通过其他行线性组合得到,它在消元过程中就会变成全零行。剩下的非零行代表了矩阵矩阵中彼此线性独立的“有效信息”数量,这个数量正是矩阵的秩。
下面给出一个具体的例子:
\[\begin{cases} x_1 + 2x_2 - x_3 = -6 \\ 2x_1 + x_2 - 3x_3 = -9 \\ -x_1 - x_2 + 2x_3 = 7 \end{cases} \]
第一步:写出增广矩阵
\[\begin{bmatrix} 1 & 2 & -1 & -6 \\ 2 & 1 & -3 & -9 \\ -1 & -1 & 2 & 7 \end{bmatrix} \]
第二步:处理第 1 列
\[\begin{bmatrix} 2 & 1 & -3 & -9 \\ 1 & 2 & -1 & -6 \\ -1 & -1 & 2 & 7 \end{bmatrix} \]
\[\begin{bmatrix} 1 & 0.5 & -1.5 & -4.5 \\ 1 & 2 & -1 & -6 \\ -1 & -1 & 2 & 7 \end{bmatrix} \]
\[\begin{bmatrix} 1 & 0.5 & -1.5 & -4.5 \\ 0 & 1.5 & 0.5 & -1.5 \\ 0 & -0.5 & 0.5 & 2.5 \end{bmatrix} \]
第三步:处理第 2 列
\[\begin{bmatrix} 1 & 0.5 & -1.5 & -4.5 \\ 0 & 1 & 1/3 & -1 \\ 0 & -0.5 & 0.5 & 2.5 \end{bmatrix} \]
\[\begin{bmatrix} 1 & 0 & -5/3 & -4 \\ 0 & 1 & 1/3 & -1 \\ 0 & 0 & 2/3 & 2 \end{bmatrix} \]
第四步:处理第 3 列
\[\begin{bmatrix} 1 & 0 & -5/3 & -4 \\ 0 & 1 & 1/3 & -1 \\ 0 & 0 & 1 & 3 \end{bmatrix} \]
\[\begin{bmatrix} 1 & 0 & 0 & 1 \\ 0 & 1 & 0 & -2 \\ 0 & 0 & 1 & 3 \end{bmatrix} \]
最终结果:
\(x_1 = 1, x_2 = -2, x_3 = 3\)。
double a[N][N]; // 增广矩阵 n * (n + 1)
double ans[N]; // 存储解
int gauss(int n) { // n: 变量个数
int r = 0; // 当前处理的行
for (int c = 0; c < n; ++c) { // 处理每一列
int t = r;
for (int i = r; i < n; ++i) // 找当前列绝对值最大的行(选主元)
if (fabs(a[i][c]) > fabs(a[t][c])) t = i;
if (fabs(a[t][c]) < EPS) continue; // 该列全为0,跳过
if (t != r) {
for (int i = 0; i <= n; i++) swap(a[t][i], a[r][i]]);
}
// 将当前行首位化为1
for (int i = n; i >= c; --i) a[r][i] /= a[r][c];
// 用当前行消去其他所有行的第 c 列
for (int i = 0; i < n; ++i) {
if (i != r && fabs(a[i][c]) > EPS) {
double f = a[i][c];
for (int j = c; j <= n; ++j)
a[i][j] -= a[r][j] * f;
}
}
r++;
}
if (r < n) {
for (int i = r; i < n; ++i) {
if (fabs(a[i][n]) > EPS) return 1; // 0 = d 形式
}
return 2; // 0 = 0 形式
}
for (int i = 0; i < n; ++i) ans[i] = a[i][n];
return 0;
}
由于有三层循环(枚举列、枚举行、消元),时间复杂度为 \(O(n^3)\)。
给定包含 \(n\) 个方程和 \(n\) 个未知数的线性方程组,要求求解该方程组。如果存在唯一解,输出各未知数的值(保留两位小数);否则输出
No Solution。
模板题。
#include <cstdio>
#include <cmath>
#include <algorithm>
using namespace std;
const int N = 105;
const double EPS = 1e-6;
// a[i][j] 存储增广矩阵,第 i 行第 j 列系数
// a[i][n] 存储等号右侧的常数 b[i]
double a[N][N];
bool solve(int n) {
int r = 0; // 当前处理的行索引
for (int c = 0; c < n; c++) { // 枚举每一列
int t = r;
// 寻找主元:在当前列及以下寻找绝对值最大的行,以保证数值稳定性
for (int i = r; i < n; i++) {
if (fabs(a[i][c]) > fabs(a[t][c])) t = i;
}
// 如果最大主元接近 0,说明该列所有元素都接近 0,此列无法消元,跳过
if (fabs(a[t][c]) < EPS) continue;
// 交换当前行与主元所在的行
if (t != r) {
for (int i = 0; i <= n; i++) swap(a[t][i], a[r][i]);
}
// 将当前行的主元系数化为 1
for (int i = n; i >= c; i--) a[r][i] /= a[r][c];
// 消去其他所有行在该列的系数(包括上方和下方的行)
for (int i = 0; i < n; i++) {
if (i != r && fabs(a[i][c]) > EPS) {
double f = a[i][c];
for (int j = c; j <= n; j++) {
a[i][j] -= f * a[r][j];
}
}
}
r++;
}
if (r < n) return -1;
else return 0;
}
int main()
{
int n; scanf("%d", &n);
for (int i = 0; i < n; i++) {
for (int j = 0; j <= n; j++) {
scanf("%lf", &a[i][j]);
}
}
if (solve(n) != 0) {
printf("No Solution\n");
} else {
// 输出解,保留两位小数
for (int i = 0; i < n; i++) {
printf("%.2f\n", a[i][n]);
}
}
return 0;
}
模板题,和上一题相比需要细分三种情况。
#include <cstdio>
#include <cmath>
#include <algorithm>
using namespace std;
const int N = 55;
const double EPS = 1e-6;
double a[N][N]; // a[i][j] 存储增广矩阵
// 求解线性方程组,返回值:
// 0: 唯一解;-1: 无解;1: 无穷多解
int solve(int n) {
int r = 0; // 当前处理的行
for (int c = 0; c < n; c++) { // 枚举每一列(变量)
int t = r;
// 选主元:在当前行及以下寻找该列绝对值最大的行
for (int i = r; i < n; i++) {
if (fabs(a[i][c]) > fabs(a[t][c])) t = i;
}
if (fabs(a[t][c]) < EPS) continue;
// 交换当前行与主元行
if (t != r) {
for (int i = 0; i <= n; i++) swap(a[t][i], a[r][i]);
}
// 归一化:将主元系数化为 1
double d = a[r][c];
for (int i = c; i <= n; i++) a[r][i] /= d;
// 消元:利用当前行消去其他所有行的第 c 列
for (int i = 0; i < n; i++) {
if (i != r && fabs(a[i][c]) > EPS) {
double f = a[i][c];
for (int j = c; j <= n; j++) {
a[i][j] -= f * a[r][j];
}
}
}
r++; // 处理完一个有效行
}
// 检查是否有解
if (r < n) {
for (int i = r; i < n; i++) {
if (fabs(a[i][n]) > EPS) return -1; // 无解
}
return 1; // 无穷多解
}
return 0;
}
int main()
{
int n; scanf("%d", &n);
for (int i = 0; i < n; i++) {
for (int j = 0; j <= n; j++) {
scanf("%lf", &a[i][j]);
}
}
int res = solve(n);
if (res == -1) printf("-1\n");
else if (res == 1) printf("0\n");
else {
// 唯一解:输出每个变量的结果
for (int i = 0; i < n; i++) {
printf("x%d=%.2f\n", i + 1, a[i][n]);
}
}
return 0;
}
给定 \(n\) 维空间球面上的 \(n+1\) 个点,要求求出该球面的球心坐标。
设球心坐标为 \(X = (x_1, x_2, \dots, x_n)\),球面上第 \(i\) 个点的坐标为 \(P_i = (p_{i,1}, p_{i,2}, \dots, p_{i,n})\)。
根据球心的定义,球心到球面上所有点的距离相等(等于半径 \(R\))。因此,对于每一个点 \(i \in [1,n+1]\),都有 \(\sum \limits_{j=1}^n (p_{i,j} - x_j)^2 = R^2\)。
这是一个包含 \(n+1\) 个方程和 \(n+1\) 个未知数(\(x_1, \dots x_n\) 以及 \(R\))的二次方程组,为了求解,需要对其进行线性化。
通过将相邻的两个方程相减(例如第 \(i+1\) 个方程减去第 \(i\) 个方程),可以消去二次项 \(x_j^2\) 和 \(R^2\)。
\[\sum_{j=1}^n (p_{i+1,j} - x_j)^2 - \sum_{j=1}^n (p_{i,j} - x_j)^2 = 0 \]
展开并简化:
\[\sum_{j=1}^n (p_{i+1,j}^2 - 2p_{i+1,j}x_j + x_j^2 - (p_{i,j}^2 - 2p_{i,j}x_j + x_j^2)) = 0 \]
\[\sum_{j=1}^n (p_{i+1,j}^2 - p_{i,j}^2 - 2x_j(p_{i+1,j} - p_{i,j})) = 0 \]
整理得到线性方程组形式:
\[\sum_{j=1}^n 2(p_{i+1,j} - p_{i,j})x_j = \sum_{j=1}^n (p_{i+1,j}^2 - p_{i,j}^2) \]
这样就得到了一个包含 \(n\) 个方程和 \(n\) 个未知数 \(x_j\) 的线性方程组 \(Ax=b\),其中 \(A_{i,j} = 2(p_{i+1,j} - p_{i,j})\),\(B_i = \sum \limits_{j=1}^n (p_{i+1,j}^2 - p_{i,j}^2)\)。
于是便可通过高斯消元法求解。
#include <cstdio>
#include <cmath>
#include <algorithm>
using namespace std;
double p[15][15], m[15][15];
int main()
{
int n; scanf("%d", &n);
// 读取 n+1 个点的坐标
for (int i = 1; i <= n + 1; i++) {
for (int j = 1; j <= n; j++) {
scanf("%lf", &p[i][j]);
}
}
// 构建线性方程组 Ax = B
// 每一行 i 是由第 i+1 个点和第 i 个点的距离公式相减得到的
for (int i = 1; i <= n; i++) {
for (int j = 1; j <= n; j++) {
m[i][j] = 2.0 * (p[i + 1][j] - p[i][j]);
m[i][n + 1] += p[i + 1][j] * p[i + 1][j] - p[i][j] * p[i][j];
}
}
// 高斯-约旦消元法
for (int i = 1; i <= n; i++) {
// 主元选择
int pe = i;
for (int j = i + 1; j <= n; j++) {
if (fabs(m[j][i]) > fabs(m[pe][i])) {
pe = j;
}
}
// 交换行
if (pe != i) {
for (int j = 1; j <= n + 1; j++) {
swap(m[i][j], m[pe][j]);
}
}
// 归一化当前行
double d = m[i][i];
for (int j = i; j <= n + 1; j++) {
m[i][j] /= d;
}
// 消去其他行的当前列
for (int j = 1; j <= n; j++) {
if (i != j) {
double f = m[j][i];
for (int k = i; k <= n + 1; k++) {
m[j][k] -= f * m[i][k];
}
}
}
}
// 输出结果,精确到小数点后 3 位
for (int i = 1; i <= n; i++) {
printf("%.3f ", m[i][n + 1]);
}
return 0;
}
有 \(n\) 个开关,初始处于开或关的状态,每个开关的操作会翻转自身以及与之关联的其他开关的状态,求这 \(n\) 个开关从初始状态到达目标状态的开关操作方案数(每个开关最多允许操作一次)。
由于状态只有“开”和“关”(可以用 \(0\) 和 \(1\) 表示),且操作的效果是“翻转”,这本质上是一个关于异或运算的线性方程组。
设 \(x_i \in \{0,1\}\) 表示是否操作第 \(i\) 个开关(\(1\) 表示操作,\(0\) 表示不操作),设 \(s_i\) 为第 \(i\) 个开关的初始状态,\(t_i\) 为目标状态。设 \(A_{j,i}=1\) 表示操作开关 \(i\) 会影响开关 \(j\) 的状态,否则 \(A_{j,i}=0\)。根据题目,操作开关 \(i\) 自身一定会改变其状态,故 \(A_{i,i}=1\)。
对于每一个开关 \(j\),其最终状态由初始状态和所有影响它的操作共同决定:
\[s_j \oplus (A_{j,1}x_1 \oplus A_{j,2}x_2 \oplus \dots \oplus A_{j,n}x_n) = t_j \]
整理得:
\[\bigoplus_{i=1}^{n} A_{j,i}x_i = s_j \oplus t_j \]
这是一个含有 \(n\) 个变量、\(n\) 个方程的异或线性方程组,其同样可以直接使用高斯消元法来求解。
如果消元后出现形如 \(0=1\) 的矛盾方程(系数全为 \(0\) 但常数项为 \(1\)),则输出 Oh,it's impossible~!!。设矩阵的秩为 \(r\),则自由元的个数为 \(n-r\)。由于每个变量只有 \(0\) 或 \(1\) 两种取值,每个自由元都可以独立选择,因此总方案数为 \(2^{n-r}\)。
时间复杂度为 \(O(n^3)\),配合 bitset 优化可达到 \(O(n^3/w)\)。
#include <cstdio>
#include <bitset>
#include <algorithm>
using namespace std;
const int N = 35;
int s[N], t[N];
bitset<N> m[N]; // m[i] 存储第 i 个方程的系数以及常数项(在第 n 位)
// 解决单组测试数据
// 问题转化为异或线性方程组 Ax = B
// A[j][i] = 1 表示开关 i 的操作会影响开关 j 的状态
// B[j] = 初始状态[j] XOR 目标状态[j]
void solve() {
int n; scanf("%d", &n);
for (int i = 0; i < n; i++) scanf("%d", &s[i]);
for (int i = 0; i < n; i++) scanf("%d", &t[i]);
for (int i = 0; i < n; i++) {
m[i].reset();
m[i][i] = 1; // 每个开关操作都会影响它自己
if (s[i] ^ t[i]) m[i][n] = 1; // 方程右侧的常数项
}
while (true) {
int i, j; scanf("%d%d", &i, &j);
if (i == 0 && j == 0) break;
// 如果操作第 i 个开关,第 j 个开关的状态也会变化
// 即 x[i] 是第 j 个方程的一个项
m[j - 1][i - 1] = 1;
}
// 高斯-约旦消元法
int r = 0;
for (int c = 0; c < n; c++) {
int p = -1;
for (int i = r; i < n; i++) {
if (m[i][c]) {
p = i;
break;
}
}
if (p != -1) {
swap(m[r], m[p]);
for (int i = 0; i < n; i++) {
if (i != r && m[i][c]) m[i] ^= m[r];
}
r++;
}
}
// 检查方程组是否相容(是否有矛盾方程 0 = 1)
for (int i = r; i < n; i++) {
if (m[i][n]) {
printf("Oh,it's impossible~!!\n");
return;
}
}
// 解的数量为 2^(自由元个数)
int ans = 1 << (n - r);
printf("%d\n", ans);
}
int main()
{
int k; scanf("%d", &k);
while (k--) solve();
return 0;
}
工厂生产 \(N\) 种零件,每种零件的生产天数 \(x_i \in [3, 9]\) 且为整数。现有 \(M\) 条记录,每条记录包含:
- 生产开始和结束的星期(SUN, MON, TUE, WED, THU, FRI, SAT)。
- 这段时间内生产的零件总数 \(k\) 以及具体的零件种类。
要求根据这些记录求出每种零件的生产天数。若无解输出
Inconsistent data.,若有多解输出Multiple solutions.,否则输出唯一的解序列。
设第 \(i\) 种零件的生产天数为 \(x_i\),对于每一条记录,假设开始星期为 \(S\),结束星期为 \(T\),零件 \(i\) 出现了 \(c_i\) 次,生产总天数(含起止当天)可以表示为 \((T - S + 1) \pmod 7\)。
注意:这里的星期需要映射到 \(0 \sim 6\),例如:MON 为 1,TUE 为 2,...,SUN 为 0。
方程可写作:
\[\sum_{i=1}^{N} c_i \times x_i \equiv (T - S + 1) \pmod 7 \]
由于零件的生产天数在 \([3, 9]\) 之间,长度正好为 7。这意味着一旦确定了 \(x_i \pmod 7\) 的值,其在 \([3, 9]\) 范围内的取值是唯一的。
在模 7 意义下进行高斯消元,因为 7 是质数,可以直接使用逆元进行除法运算。
#include <cstdio>
#include <algorithm>
using namespace std;
const int N = 305;
const int INV[] = {0, 1, 4, 5, 2, 3, 6}; // 模 7 逆元表
int n, m, a[N][N];
char s[4], t[4];
int get(char s[]) {
if (s[0] == 'M') return 1;
if (s[0] == 'T' && s[1] == 'U') return 2;
if (s[0] == 'W') return 3;
if (s[0] == 'T' && s[1] == 'H') return 4;
if (s[0] == 'F') return 5;
if (s[0] == 'S' && s[1] == 'A') return 6;
if (s[0] == 'S' && s[1] == 'U') return 0;
return -1;
}
void solve() {
// 初始化矩阵
for (int i = 0; i < N; i++)
for (int j = 0; j < N; j++)
a[i][j] = 0;
for (int i = 0; i < m; i++) {
int k;
scanf("%d%s%s", &k, s, t);
int tg = (get(t) - get(s) + 1 + 7) % 7;
for (int j = 0; j < k; j++) {
int tp; scanf("%d", &tp);
a[i][tp - 1] = (a[i][tp - 1] + 1) % 7;
}
a[i][n] = tg;
}
// 高斯消元
int r = 0, c = 0;
while (r < m && c < n) {
int p = r;
for (int i = r; i < m; i++) {
if (a[i][c] != 0) {
p = i;
break;
}
}
if (a[p][c] == 0) {
c++;
continue;
}
for (int j = 0; j <= n; j++) swap(a[r][j], a[p][j]);
for (int i = 0; i < m; i++) {
if (i != r && a[i][c] != 0) {
int f = a[i][c] * INV[a[r][c]] % 7;
for (int j = c; j <= n; j++) {
a[i][j] = (a[i][j] - f * a[r][j] % 7 + 7) % 7;
}
}
}
r++; c++;
}
// 检查一致性
for (int i = r; i < m; i++) {
if (a[i][n] != 0) {
printf("Inconsistent data.\n");
return;
}
}
// 检查多解
if (r < n) {
printf("Multiple solutions.\n");
return;
}
// 提取唯一解并映射到 [3, 9]
for (int i = 0; i < n; i++) {
int res = a[i][n] * INV[a[i][i]] % 7;
if (res < 3) res += 7;
printf("%d%c", res, i == n - 1 ? '\n' : ' ');
}
}
int main()
{
while (true) {
scanf("%d%d", &n, &m);
if (n == 0 && m == 0) break;
solve();
}
return 0;
}
此内容由惯性聚合(RSS阅读器)自动聚合整理,仅供阅读参考。 原文来自 — 版权归原作者所有。