高斯消元学习笔记

$n$ 元一次方程组

使用场景

高斯消元一般是用来求解线性方程组的方法。

线性同余方程组大家都见过，相比线性方程组应该并不陌生。

至于线性方程组，是线性代数中的一个很小的部分。是指类似这样的一个方程组：

$$\{\begin{array}{l}{a}_{1,1}{x}_1+a_{1,2}{x}_2+\cdots +a_{1,n}{x}_n=b_1\\ a_{2,1}{x}_1+a_{2,2}{x}_2+\cdots +a_{2,n}{x}_n=b_2\\ \cdots \\ a_{n,1}{x}_1+a_{n,2}{x}_2+\cdots +a_{n,n}{x}_n=b_n\end{array}$$

也就是有 $n$ 个未知数（{x1,x2,⋯ ,xn}\{x_1,x_2,\cdots,x_n\}{x1​,x2​,⋯,xn​}）和 $n$ 个方程（{L1,L2,⋯ ,Ln}\{L_1,L_2,\cdots,L_n\}{L1​,L2​,⋯,Ln​}）。方程的样子就是每一个未知数带上一个系数，加起来得到另一个数。

也有可能是 $n$ 个未知数和 $m$ 个方程（$n\ne m$），但是 $m<n$ 的话就可能是不定方程，$m>n$ 的话也没啥必要。所以就只需要考虑 $n=m$ 即可。

---

算法思路

那么这东西到底是怎么做的呢？

首先，考虑将每一个未知数的系数和结果都拎出来，因为我们已经知道 $x_1,x_2,\cdots ,x_n$ 是有顺序排列的未知数没有用，所以只需要将这个方程的“灵魂”弄出来就行了。

考虑将每一个方程 $L_i$ 的所有未知数的系数和 $L_i$ 的右侧排在一起。$a_{i,1}$ 排在第一个，$a_{i,2}$ 排在第二个……$a_{i,n}$ 排在倒数第二个，$b_i$ 排在最后一个。

最终会形成一个矩阵。

---

例如

$$\{\begin{array}{l}2x+y-z=8\\ -3x-y+2z=-11\\ -2x+y+2z=-3\end{array}$$

这个线性方程组变成的矩阵是：

$$\left[\begin{array}{cccc}2& 1& -1& 8\\ -3& -1& 2& -11\\ -2& 1& 2& -3\end{array}\right]$$

于是我们就做完了高斯消元的第一步。

---

模型已经建立好了。

高斯消元的目的就是：将这个矩阵前面 $n\times n$ 的系数构成的矩阵，变成一个上三角矩阵。也就是副对角线的下面全是 $0$，上面才会可能有值。

即

$$\left[\begin{array}{cccc}2& 1& -1& 8\\ -3& -1& 2& -11\\ -2& 1& 2& -3\end{array}\right]$$

会变成

$$\left[\begin{array}{cccc}?& ?& ?& 8\\ 0& ?& ?& -11\\ 0& 0& ?& -3\end{array}\right]$$

其中问号代表的是不确定的部分，因为显然消掉副对角线下面的系数方程会发生变化。

所以，高斯消元的“消元”指的就是把这些东西消掉。

---

这下高斯消元的意图就很明显了吧——这样就可以倒着求出所有的解。

显然这样就能较为简单地求解。

因为这个时候 $L_n$ 是一个十分简单的形如 $ax=b$ 的方程，可以很容易求解 $x_n$ 的值。

然后将其代入上面一个式子，即又可以得到 $L_{n-1}$ 也是一个 $ax=b$ 的方程，又可以得出 $x_{n-1}$ ……

就这样一直代入，我们就可以得出所有未知数的解。

---

那么它到底是怎么消元的呢？？实际上运用的还是我们初中数学学的消元法。

还记得消元法怎么弄的吗？？？就是其中一个方程乘上一个系数，加上另外一个方程乘上一个系数，得到另外一个方程。然后每一个未知数的系数发生改变，右边的值也发生改变。特别地，其中还有一个未知数的系数变成了 $0$，这就完成了消元的操作。

首先我们考虑分步处理：先处理第一列的 $0$，再处理第二列的 $0$……

这个处理过程，转换到实际上的方程就变成：先将 $L_2,\cdots ,L_n$ 的 $x_1$ 系数变成 $0$（即将 $x_1$ 消掉），再将 $L_3,\cdots ,L_n$ 的 $x_2$ 消掉……

不妨再次使用文章中反复提到的

$$\{\begin{array}{l}2x+y-z=8\\ -3x-y+2z=-11\\ -2x+y+2z=-3\end{array}$$

这个线性方程组来演示消元的过程。

---

首先，考虑处理第一列。

我们可以：

$$\{\begin{array}{l}{L}_2+\frac{3}{2}{L}_1\to L_2\\ L_3+L_1\to L_3\end{array}$$

这样就可以消掉了。

这个时候方程就变为：

$$\{\begin{array}{l}2x+y-z=8\\ \frac{1}{2}y+\frac{1}{2}z=1\\ 2y+z=5\end{array}$$

于是我们顺利地消掉了第一列的某一些系数。在这个时候要消元的方程都是使用了 $L_1$ 来消元。

---

考虑处理第二列的 $0$。

这个时候我们也只需要对 $L_3$ 下手脚了，因为 $n=3$。

如果 $L_3$ 加上 $?\times L_1$ 的话 $x$ 又会回来，于是这次必须要使用 $L_2$ 来消元。

于是我们可以：

$$L_3+(-4)L_2\to L_3$$

这样的方程就变为：

$$\{\begin{array}{l}2x+y-z=8\\ \frac{1}{2}y+\frac{1}{2}z=1\\ -z=1\end{array}$$

至此，这个例子的消元完毕了。然后就是按部就班地得出 $x,y,z$ 的结果：

$$\{\begin{array}{l}x=2\\ y=3\\ z=-1\end{array}$$

然后再代入回原来的方程来计算，嗯，正确的。

---

在这里的时候我们会在下面选未知数的系数绝对值最大的方程与其交换，因为这个时候高斯消元可能会出现浮点数，会出现精度问题，所以选一个最大的被除数以减小数值计算的精度误差。

所以，我们的举例使用的做法可能不太对，不管怎样，都要保证消元使用的未知数的系数的绝对值最大。

---

你有没有发现一个规律，对于第 $i$ 列要消去 $x_i$ 的方程 $L_{i+1},\cdots ,L_n$，必须要使用方程 $L_i$ 来进行消元。

而且很容易发现，当 $L_x$ 要使用 $L_i$ 这个方程来消元的时候，都是 $-\frac{a_{x,i}}{a_{i,i}}{L}_i+L_x\to L_x$。

---

总结一下，高斯消元得到答案的步骤主要有以下三个步骤：

• 首先将所有方程里面的系数和结果的值都拎出来，组成一个矩阵。

• 一列一列地考虑消元，套用上面的公式即可，不过还要注意顺带更新其他未知数的系数。

• 从最后一行开始，求出 $x_n,x_{n-1},\cdots ,x_1$ 的值。

很容易发现，这种 $n$ 个未知数 $n$ 个方程的线性方程组一般来说仅会有一组解。

无解判断

显然这个时候有可能会出现无解的情况或者是无穷多组解的情况。

无解情况：高斯消元完了之后，从最后使用来消元的方程到最后一个方程的 $a$ 理应都全部是 $0$。如果这个时候有一个 $b\ne 0$ 不久废了吗！

无穷解情况：高斯消元完了之后，发现还有一些方程没有被涉及到，而且不是无解情况，显然这些方程的 $a$ 都应该是 $0$，而就会有无数多组解（我取那个数都可以，反正有 $0$，答案就一定是 $0$）。

代码

以 P3389 为例子。

#include <bits/stdc++.h>
using namespace std;

const int N = 110;          // 最大方程数和变量数
const double eps = 1e-10;    // 浮点数精度阈值

int n;                       // 方程数和变量数
double a[N][N], b[N];        // a为系数矩阵，b为常数项数组

// 高斯-约旦消元法解线性方程组，转化为对角矩阵直接求解
void gauss() {
    int l = 1;               // 当前主元所在的行
    for (int i = 1; i <= n; i++) {  // 枚举每一列
        // 寻找第i列中绝对值最大的行，进行行交换（选主元）
        for (int j = l; j <= n; j++) {
            if (fabs(a[j][i]) > fabs(a[l][i])) {
                // 交换第l行和第j行
                for (int k = i; k <= n; k++)
                    swap(a[l][k], a[j][k]);
                swap(b[l], b[j]);
            }
        }
        // 若主元接近0，跳过该列（视为自由变量，但此处可能导致无解或无穷解）
        if (fabs(a[l][i]) < eps)
            continue;
        // 用主元行消去其他行的第i列元素
        for (int j = 1; j <= n; j++) {
            if (j != l && fabs(a[j][i]) > eps) {
                double d = a[j][i] / a[l][i]; // 计算消元系数
                // 消元：行j -= 行l * d
                for (int k = i; k <= n; k++)
                    a[j][k] -= a[l][k] * d;
                b[j] -= b[l] * d;
            }
        }
        l++; // 主元行处理完毕，移至下一行
    }
    // 检查无解情况：存在0x=非零的方程
    for (int i = l; i <= n; i++) {
        if (fabs(b[i]) > eps) {
            cout << "No Solution\n";
            return;
        }
    }
    // 检查是否有足够的主元（唯一解条件）
    if (l <= n) { // 主元数量不足，可能存在自由变量，无穷解视为无解
        cout << "No Solution\n";
    } else {
        // 输出唯一解（此时矩阵已转为对角形，直接计算）
        for (int i = 1; i <= n; i++)
            printf("%.2lf\n", b[i] / a[i][i]);
    }
}

int main() {
    cin >> n;
    for (int i = 1; i <= n; i++) {
        for (int j = 1; j <= n; j++)
            cin >> a[i][j];
        cin >> b[i];
    }
    gauss();
    return 0;
}

复杂度 $O(n^3)$。