高斯消元

高斯消元 \(O(n^3)\)

对于一个 \(n*(n+1)\) 的矩阵, 有 :

\[\begin{aligned} a_{11}x_1+a_{12}x_2+&...+a_{1n}x_n=b_1 \\ a_{21}x_1+a_{22}x_2+&...+a_{2n}x_n=b_2 \\ . \\ . \\ . \\ a_{n1}x_1+a_{n2}x_2+&...+a_{nn}x_n=b_n \end{aligned} \]

可以解出 \(\{x_1,x_2,x_3,...,x_n\}\).

算法步骤如下 :

1. 找主元 : 找到当前列中系数绝对值最大的一行.
2. 归一化 : 将改行交换到最上面 (未确定的一行).
3. 消 : 用该行将下面所有行的当且列的值消为 \(0\).

做完上面的步骤, 矩阵已经被化为了上三角型.

接下来判断解的情况 :

1. 出现了 \(0\ !=0\), 说明无解.
2. 未出现 (1), 但出现了 \(0=0\), 说明有无穷多组解.
3. 未出现 (1), (2), 说明解唯一.

下面是代码模板 :

int gauss() {
    int c, r;
    for (c = 0, r = 0; c < n; ++c) {
        int t = r;
        for (int i = r + 1; i < n; ++i)
            if (fabs(a[i][c]) > fabs(a[t][c]))
                t = i;
        if (fabs(a[t][c]) < eps) continue;
        for (int i = c; i <= n; ++i) swap(a[t][i], a[r][i]); // 交换
        for (int i = n; i >= c; --i) a[r][i] /= a[r][c];  // 归一化
        for (int i = r + 1; i < n; ++i) // 消
            if (fabs(a[i][c]) > eps)
                for (int j = n; j >= c; --j) a[i][j] -= a[r][j] * a[i][c];
        ++r;
    }
    if (r < n) {
        for (int i = r; i < n; ++i)
            if (fabs(a[i][n]) > eps) return 2; // 0 != 0, 无解
        return 1; // 无穷多解
    } 
   	// 唯一解, 接下来确认解的具体值, 即把矩形化为对角线型
    for (int i = n - 1; i >= 0; --i)
        for (int j = i + 1; j < n; ++j)
            a[i][n] -= a[j][n] * a[i][j];
    return 0;
}