zoukankan      html  css  js  c++  java
  • 「笔记」高斯消元

    消元法

    先来看一下百度百科的定义:

    消元法是指将许多关系式中的若干个元素通过有限次地变换,消去其中的某些元素,从而使问题获得解决的一种解题方法。

    可能不好懂。

    回想一下小学数学中解二元一次方程的方法

    比如下面这个二元一次方程:

    [egin{cases} x + y = 10\ x - y = 6 end{cases} ]

    通过两式相加,可以得到

    [2x = 16 ]

    进而解得

    [x=8 ]

    (x) 代回就可以得到

    [y = 2 ]

    发现,在解这个二元一次方程时,通过式子加减,我们把其中一个元 (y) 消掉,只剩下了 (x),变成了我们可以直接得出解的形式。

    这个解题方法就是消元法

    高斯消元其实就是在模拟这个消元过程。

    只不过元可能更多,式子可能更复杂,但好在这是计算机考虑的问题,与我们没有什么关系。

    高斯消元

    高斯-约旦消元法

    与后面高斯消元法相比没有回带过程,实现起来也更加简单,精度方面也更加优秀

    高斯消元的一般步骤是:

    • 列出要进行消元的矩阵
    • 枚举每一列
    • 设枚举到第 (i) 列,找到第 (i) 列中系数最大的方程并将其交换到第 (i)
    • 同过确定系数,加减消元等方法,将其他方程中的这个元的系数化为 (0)
    • 求出每个未知数的解

    进行高斯消元后是一个三角矩阵,而高斯-约旦消元法的结果是一个对角线矩阵

    举个例子:

    求解下面这个方程组

    规定从上到下三个式子分别是 (A), (B), (C)

    [egin{cases} 3x + 2y + z = 10 \ 5x + y +6z = 25 \ 2x + 3y + 4z = 20 end{cases} ]

    把它转化成矩阵形式

    [egin{bmatrix} 3 & 2 & 1 & mid & 10 \ 5 & 1 & 6 & mid & 25 \ 2 & 3 & 4 & mid & 20 end{bmatrix} ]

    左边是每个未知数的常数项,右边是该式结果,中间用一条竖线隔开

    我们从第一列开始消元

    找到第一列最大的一个,与第一行交换

    [egin{bmatrix} 5 & 1 & 6 & mid & 25 \ 3 & 2 & 1 & mid & 10 \ 2 & 3 & 4 & mid & 20 end{bmatrix} ]

    我们要将第一列的其他两项化为 (0)

    (B = B - frac{3}{5}A)

    (C = C - frac{2}{5}A)

    然后变成这样:

    [egin{bmatrix} 5 & 1 & 6 & mid & 25 \ 0 & frac{7}{5} & -frac{13}{5} & mid & -5 \ 0 & frac{12}{5} & frac{2}{5} & mid & 5 end{bmatrix} ]

    然后在看第二列,把第二列系数最大的那个方程换到第二列。

    [egin{bmatrix} 5 & 1 & 6 & mid & 25 \ 0 & frac{12}{5} & frac{2}{5} & mid & 5 \ 0 & frac{7}{5} & -frac{13}{5} & mid & -5 end{bmatrix} ]

    继续进行消元操作,不过这次我们要把方程 (A) 的第二列的元也消掉

    (A = A - frac{5}{12}B)

    (C = C - frac{7}{12}B)

    矩阵变为:

    [egin{bmatrix} 5 & 0 & frac{35}{6} & mid & frac{275}{12} \ 0 & frac{12}{5} & frac{2}{5} & mid & 5 \ 0 & 0 & frac{17}{6} & mid & frac{95}{12} end{bmatrix} ]

    继续对第三列进行相同操作,矩阵变为:

    [egin{bmatrix} 5 & 0 & 0 & mid & frac{225}{34} \ 0 & frac{12}{5} & 0 & mid & 5 \ 0 & 0 & frac{17}{6} & mid & frac{66}{17} end{bmatrix} ]

    到此消元结束,发现每个方程只剩下一个变量,直接解出来输出结果即可。

    如果不懂的话可以看一下代码实现。

    /*
    Work by: Suzt_ilymics
    Problem: P3389 【模板】高斯消元法
    Knowledge: 高斯-约旦消元法 
    Time: O(n^3)
    */
    #include<iostream>
    #include<cstdio>
    #include<cstring>
    #include<algorithm>
    #include<cmath>
    #define LL long long
    #define orz cout<<"lkp AK IOI!"<<endl
    
    using namespace std;
    const int MAXN = 105;
    const int INF = 1e9+7;
    const int mod = 1e9+7;
    
    int n;
    double a[MAXN][MAXN];
    
    int read(){
        int s = 0, f = 0;
        char ch = getchar();
        while(!isdigit(ch))  f |= (ch == '-'), ch = getchar();
        while(isdigit(ch)) s = (s << 1) + (s << 3) + ch - '0' , ch = getchar();
        return f ? -s : s;
    }
    
    int main()
    {
        n = read();
        for(int i = 1; i <= n; ++i) 
            for(int j = 1; j <= n + 1; ++j) 
                scanf("%lf", &a[i][j]); // 读入方程组 
        
        for(int i = 1; i <= n; ++i) {
            int Max = i; 
            for(int j = i + 1; j <= n; ++j) if(fabs(a[i][j]) > fabs(a[i][j])) Max = j; // 选出该列最大系数 
            for(int j = 1;j <= n + 1; ++j) swap(a[i][j], a[Max][j]);  // 交换这两行 
            if(!a[i][i]) { // 最大值等于 0 则说明该列都为 0,肯定无解 
                puts("No Solution");
                return 0;
            }
            for(int j = 1; j <= n; ++j) { // 对每一行进行加减消元 
                if(j != i) {
                    double tmp = a[j][i] / a[i][i];
                    for(int k = i + 1; k <= n + 1; ++k) a[j][k] -= a[i][k] * tmp;
                }
            }
        }
        for(int i = 1; i <= n; ++i) printf("%.2lf
    ", a[i][n + 1] / a[i][i]);
        return 0;
    }
    

    高斯消元法

    没写过。

    实现起来与高斯-约旦类似,添了一个回带的过程

    在这里粘一个 皎月半撒花 大佬的板子

    #include<cstdio>
    #include<iostream>
    #include<algorithm>
    #include<cmath>
    using namespace std;
    double map[111][111];
    double ans[111];
    double eps=1e-7;
    int main(){
        int n;
        cin>>n;
        for(int i=1;i<=n;i++)
            for(int j=1;j<=n+1;j++)
                scanf("%lf",&map[i][j]);
        for(int i=1;i<=n;i++){
            int r=i;
            for(int j=i+1;j<=n;j++)
                if(fabs(map[r][i])<fabs(map[j][i]))
                    r=j;//find_the_biggest_number_of_the_first_column(at present) 
            if(fabs(map[r][i])<eps){
                printf("No Solution");
                return 0;
            }
            if(i!=r)swap(map[i],map[r]);//对换一行或一列,属于找最大当前系数的其中一步。(这样就可以只处理当前行的系数啦!) 
            double div=map[i][i];
            for(int j=i;j<=n+1;j++)
                map[i][j]/=div;
            for(int j=i+1;j<=n;j++){
                div=map[j][i];
                for(int k=i;k<=n+1;k++)
                    map[j][k]-=map[i][k]*div;
            }
        }
        ans[n]=map[n][n+1];
        for(int i=n-1;i>=1;i--){
            ans[i]=map[i][n+1];
            for(int j=i+1;j<=n;j++)
                ans[i]-=(map[i][j]*ans[j]);
        }//回带操作
        for(int i=1;i<=n;i++)
            printf("%.2lf
    ",ans[i]);
    }
    
    
  • 相关阅读:
    shell log
    Python:列出列表中所有元素的组合可能
    scrapy 停止爬虫
    shell split log by data
    mitmproxy 配置
    插件reres的使用,替换网站的js文件
    解决小米Note adb调试无法发现设备
    md5 计算文件一致性
    【Frida Hook 学习记录】Frida Hook Android 常用方法
    监控神器普罗米修斯Prometheus安装配置
  • 原文地址:https://www.cnblogs.com/Silymtics/p/14797630.html
Copyright © 2011-2022 走看看