Skip to content

高斯消元法可以求解线性方程组。所谓线性方程组,其实就是由 MN 元一次方程构成的。线性方程组的所有系数可以表示为一个 MN 列的系数矩阵,如果再加上等号右侧的常数列,则可以表示为一个 MN+1 列的增广矩阵

例如,有如下线性方程组:

{x+2yz=62x+1y3z=9x+y+2z=7

可以用矩阵形式表示为:

[121213112][xyz]=[697]

其中,左侧为系数矩阵,右侧为常数列。

求解该类方程的步骤可以概括为对增广矩阵进行初等行变换。而初等行变换有以下三种:

  1. 某行乘以一个非零常数。
  2. 把其中一行的若干倍加到另一行。
  3. 交换两行的位置。

同理还有初等列变换。

上例的增广矩阵可以进行初等行变换如下:

[121621391127][213912161127][10.51.54.512161127][10.51.54.501.50.51.500.50.52.5][10.51.54.50113100.50.52.5][10.51.54.50113100232]

最后得到的矩阵被称为阶梯形矩阵,它的系数矩阵被称为上三角矩阵,这个矩阵对应的线性方程组为:

{x+0.5y1.5z=4.5y+13z=123z=2

这样我们就求得了最后一个未知量的值,从下往上依次带入方程组,即可求得其他未知量的值。

其实该矩阵也可进一步简化,得到更简洁的形式。

[10.51.54.50113100232][100101020013][100010001][xyz]=[123]

最后得到的矩阵称为简化阶梯形矩阵,它的系数矩阵称为对角矩阵

这种方法即是高斯消元法

考虑到线性方程组的解的有三种情况:

  1. 若存在系数全为零、常数不为零的行,则方程组无解
  2. 若系数不全为零的行恰好有 N 个,则方程组唯一解
  3. 若系数不全为零的行有 k<N 个,则方程组无穷多组解
cpp
#include <bits/stdc++.h>
using namespace std;

const int N = 110;
const double eps = 1e-8;

int n;
double a[N][N];

int gauss() {
    int c, r;
    for (c = 0, r = 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;

        // 交换主元到顶端
        for (int i = c; i <= n; i++) swap(a[t][i], a[r][i]); 
        // 主元变为 1
        for (int i = n; i >= c; i--) a[r][i] /= a[r][c];
        // 主元列下的元素消为 0
        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;
        }
        
        // 无穷多组解
        return 1; 
    }

    for (int i = n - 1; i >= 0; i--) {
         for (int j = i + 1; j < n; j++) {
            a[i][n] -= a[i][j] * a[j][n];
         }
    }
    
    // 唯一解
    return 0; 
}


int main() {
    scanf("%d", &n);
    for (int i = 0; i < n; i ++ ) {
        for (int j = 0; j < n + 1; j ++ ) {
            scanf("%lf", &a[i][j]);
        }
    }

    int t = gauss();
    if (t == 2) puts("No solution");
    else if (t == 1) puts("Infinite group solutions");
    else {
        for (int i = 0; i < n; i ++ ) printf("%.2lf\n", a[i][n]);
    }

    return 0;
}