小指针

V1

2022/03/28阅读:21主题:凝夜紫

高斯消元解线性方程组

封面原图
封面原图

高斯消元法

高斯消元法是数学上线性代数中的一种算法,可以把矩阵转换成行阶梯型矩阵

高斯消元可以用来找到形如下列的多元线性方程组的解:

对于上述方程组,我们的目标是求出所有的未知数x

初等行列变换

我们如果将所有的系数拿出来就会形成一个系数矩阵,高斯消元法对矩阵做初等行列变换。

初等行列变换以三种基本变换形式构成:

  1. 把某一行乘以一个非0的数;
  2. 交换某两行;
  3. 把某行的若干倍加到另外一行;

这三种变换都是等价变换。既,变换后不会对原方程的解产生任何影响

使用初等行列变换后,最理想的情况是,我们把矩阵转换成了阶梯型矩阵:

如何求解

当变成如上形式之后,我们发现最后一项只有 一个未知数,那么就可以直接求出来啦。

再看倒数第二行,只有两个未知数,既 ,而其中的 已经被求出来了,所以在倒数第二行中我们又求出了 。以此类推,当计算完第一行之后,我们就求出了所有的未知数x

从上面的求解过程中来看,显然会有三种结果:

  1. 转换后是完美阶梯形,存在唯一解;
  2. 转换后求解过程中出现了左右明显矛盾,如:左面为0,右面不为0。此时无解。
  3. 转换后求解过程中,某个等式出现0=0的情况,说明这行的等式可以取任意数,那么也就是最终方程组会有无穷多的解。

举个例子

当前有如下线性方程组:

转换成去掉系数的矩阵:

高斯消元法
  1. 枚举每一列。找到当前列(目前当前列是第一列),绝对值最大的那一行;

  2. 将这一行换到最上面。得到:

  3. 将这行的第一个数变成1。也就是2*(1/2)=1,那么该行所有的数都要乘上1/2。所以,我们就得到了:

  4. 将下面所有行的当前列变成0。那么对于第二行,它的当前列是1,我们只要用它减去第一行,就能让它变成0,于是,计算过后变成:

  5. 对于第三行,我们只要让它加上第一行就能把它的当前列变成0。于是,计算后变成:

  6. 如此,第一轮循环就完成了,并且我们固定了第一列,剩下的方程组变成了:

  7. 继续重复步骤1~5。目前两列中绝对值最大的是当前的第一列,将它的当前行变成1,既,1.5/1.5=1,那么该行的所有值都要除以1.5,计算后得到:

  8. 再把后面所有行中当前列变成0,那么只需要让下一行加上当前行的一半:-0.5+1*1/2=0,那么下一行所有的数都如此做,计算后得到:

  9. 至此,已经把所有的行都完成了转换,最后得到梯形矩阵:


    对应方程组是:

  10. 剩下的就是从下到上,依次求出所有的未知数x即可。对于最后一个方程,两面同时乘以 ,得到

  11. 对于倒数第二行,只需要让它减去 最后一行乘以 ,计算后得到:

  12. 对于第一行,先消去第二个数0.5,我们用它减去 第二行乘以0.5,得到:

  13. 对于第一行,再消去-1.5,只要让它加上 第三行乘以1.5,最后得到:

Code

输入格式
第一行包含整数 n。
接下来 n 行,每行包含 n+1 个实数,表示一个方程的 n 个系数以及等号右侧的常数。

输入样例:
3
1.00 2.00 -1.00 -6.00
2.00 1.00 -3.00 -9.00
-1.00 -1.00 2.00 7.00

输出样例:
1.00
-2.00
3.00

#include <iostream>
#include <algorithm>
#include <cmath>

using namespace std;

const int N = 110;
// 浮点数存储有误差,不能直接判0,需要判断是否小于一个极小数。
const double eps = 1e-8;
int n;
double a[N][N]; // 存储矩阵

int gauss() {
    int c, r; // c表示列,r表示行
    for (c = 0, r = 0; c < n; c ++ ) {
        int t = r;
        //找到绝对值最大的那一行。步骤1;
        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;
        // 把这一行 换到 未处理完成行 中 的 第一行。步骤2;
        for (int i = c; i <= n; i ++ ) swap(a[t][i], a[r][i]);
        // 把这行的第一个数变成1,既让该行所有的数都除以第一个数
        // 第一个数就变成了1。这里需要注意必须从后向前算,最后一个更新第一个数。
        for (int i = n; i >= c; i -- ) a[r][i] /= a[r][c];
        // 用当前行将下面所有的列消成0;步骤4,5;
        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 == 0) {
        for (int i = 0; i < n; i ++ ) {
            if (fabs(a[i][n]) < eps) a[i][n] = 0;  // 去掉输出 -0.00 的情况
            printf("%.2lf\n", a[i][n]);
        }
    }
    else if (t == 1puts("Infinite group solutions");
    else puts("No solution");
    return 0;
}

分类:

数学

标签:

数学编程

作者介绍

小指针
V1