• [未完成][知识点]高斯消元法及其扩展


    1、前言

      中学甚至是小学,了解方程的时候,我们一定是学习过高斯消元法的,或许当时只是不是这种称呼罢了。而这个知识点,在信息学上依旧可以运用,甚至有着更多的推广的功能。

    2、方程转矩阵

      我们先从最简单的求解三元一次方程入手。现存在一个方程组,我们将其写成一个矩阵:

      这个矩阵我们称之为增广矩阵(Augmented Matrix)。左侧部分为系数矩阵,最后一列是等号右边的常数列。为什么要把他转换为矩阵,因为这更加符合信息学的制式,我们直接拿二维数组存下来就行了。三元一次方程组还是很好解的,我们借用下当年的解方程方法就可以轻松地理解这些过程。

      我们固然知道,解三元一次方程组,需要进行两次消元。推广到一般的,对于n元一次,需要进行n-1次消元,我们的步骤就是从第1行起一行行处理,直至n-1行,使最终对于第i行,a[i][i]!=0 && a[j][i]==0 (j>i)。为了提高数值稳定性,每次消元之前,对于第i行,我们需要从第i+1行至第n行选取一行,设其为第x行,存在绝对值最大的a[x][i],然后交换第x行和第i行。不要问我为什么,我也不清楚提高数值稳定性具体是个什么概念,只可意会,不可言传。

    3、加减消元

      消元的核心步骤。依旧是上面那个例子。从第1行起处理,第一步我们找到了第2行满足所述条件,将其与第1行交换,得:

      然后{-3,-1,2,-11}就成为当前行了。根据当前行,依次将后面所有的行的第1位通过加减进行消除。例如第2行每个元素*(-3/2),然后与第1行相加得到{0,1/3,1/3,2/3}。

      推广到一般情况,即用第i行消去第j行的第i列,那么第j行的第k列有a[j][k]-=a[j][i]/a[i][i],其中j∈(i+1,n),k∈(i,n+1)。

      在一步步的消除之后,对于第n行,必定仅存在一个系数矩阵元素。例如还是上述例子,消元结束后,最后一列元素为{0,0,1/5,-1/5},即z/5=-1/5 => z=-1。而将z代入上一行,可以得到y=3;将y,z的值代入上一行,可以得到x=2。这个步骤我们称之为回代。这样n个未知数就都可以一步步得出来了。

    代码:

    ----------------------------------------------------------------------------------------------------

    #include<cstdio>
    #define MAXN 1005

    double abs(double o) { return o<0?-o:o; }

    void swap(double &a,double &b) { double t=a; a=b,b=t; }

    double map[MAXN][MAXN],ans[MAXN];
    int n;

    void Guass()
    {
      for (int i=1;i<=n;i++)
      {
        int o=i;
        for (int j=i+1;j<=n;j++) o=abs(map[j][i])>abs(map[o][i])?j:o;
        if (o!=i) for (int j=1;j<=n+1;j++) swap(map[i][j],map[o][j]);
        for (int j=i+1;j<=n;j++)
        {
          double x=map[j][i]/map[i][i];
          for (int k=i;k<=n+1;k++) map[j][k]-=x*map[i][k];
        }
      }
      for (int i=n;i>=1;i--)
      {
        for (int j=i+1;j<=n;j++) map[i][n+1]-=map[j][n+1]*map[i][j];
        map[i][n+1]/=map[i][i];
      }
    }

    int main()
    {
      freopen("Gauss.in","r",stdin);
      freopen("Gauss.out","w",stdout);
      scanf("%d",&n);
      for (int i=1;i<=n;i++)
        for (int j=1;j<=n+1;j++) scanf("%lf",&map[i][j]);
      Guass();
      for (int i=1;i<=n;i++) printf("%lf ",map[i][n+1]);
      return 0;
    }

    ----------------------------------------------------------------------------------------------------

    4、异或方程组

    上面所解的方程是非常常规的,高斯消元法的用途不局限于这么简单的方程。一个新内容——异或方程组,他和一般的n元一次方程组的差别在于,对于方程左侧,一般方程为元素之间相加,而异或方程组为互相异或,即形如:(a1x1)^(a2x2)^...^(anxn)=a0。因为形式的基本一致,故运算过程也没有大幅度的变化。

  • 相关阅读:
    【PAT甲级】1079 Total Sales of Supply Chain (25 分)
    CQOI2018 Day1 社交网络
    codeforces 707E Garlands (离线、二维树状数组)
    NOI2018 Day1 归程(Kruskal重构树)
    NOI2018 Day2 屠龙勇士(扩展孙子定理+multiset)
    知识点:二叉(重量)平衡树——替罪羊树
    BZOJ3065 带插入区间K小值
    知识点:斜率优化DP
    知识点:FFT详解
    博客园test(搭博客用)
  • 原文地址:https://www.cnblogs.com/jinkun113/p/4830769.html
Copyright © 2020-2023  润新知