高斯消元详解

引入：

　　dalao若赶时间，可直接阅读正文。

　　从学数学不就后，我们就开始接触方程，这东西我们再熟悉不过了。

　　这儿我们只讨论方程未知数为一次的情况。

　　一元一次方程：

　　　　例如：3x+6 = 18，显然x = 4,毫无技巧可言。

　　两元一次方程：

　　　　例如：3x+2y=36  ①,

　　　　　　　这方程看起来方程有无数组解。

　　　　　　　而想要解出某一准确解，则需要另一方程组。

　　　　　想一下：若另一方程组是带①式中某一变量的一元一次方程你会咋么解？

　　　　　　　　　若是二元一次方程呢？

　　　　这里给出 5x+y=39 为 ② 式。

　　　　小学老师曰：“将未知量系数变为相同/相反，上下相加或相减消去一变量，然后就是解简单的一元一次方程了,所以解释什么呢？”。

　　解这个方程的过程就是一个消元的过程。

　　　消元：顾名思义就是消去“元”，也就是消去方程组中的未知量。

 

正文：

　　接下来进入正题，解一般线性方程组。

偷偷告诉你：方程组中含有几个未知量，则至少需要几个方程式才能使得方程组有解（注意是至少）。

　　而解方程的人，也是今天的主角---高斯。

 

人才
高斯7岁那年开始上学。10岁的时候，他进入了学习数学的班级，这是一个首次创办的班，孩子们在这之前都没有听说过算术这么一门课程。数学教师是布特纳，他对高斯的成长也起了一定作用。
一天，老师布置了一道题，1+2+3······这样从1一直加到100等于多少。
高斯很快就算出了答案，起初高斯的老师布特纳并不相信高斯算出了正确答案：“你一定是算错了，回去再算算。
”高斯说出答案就是5050，高斯是这样算的1+100=101，2+99=101······1加到100有50组这样的数，所以50X101=5050。
布特纳对他刮目相看。他特意从汉堡买了最好的算术书送给高斯，说：“你已经超过了我，我没有什么东西可以教你了。
”接着，高斯与布特纳的助手巴特尔斯建立了真诚的友谊，直到巴特尔斯逝世。他们一起学习，互相帮助，高斯由此开始了真正的数学研究。

成就
　数学
　　欧几里得已经指出，正三边形、正四边形、正五边形、正十五边形和边数是上述边数两倍的正多边形的几何作图是能够用圆规和直尺实现的，但从那时起关於这个问题的研究没有多大进展。
高斯在数论的基础上提出了判断一给定边数的正多边形是否可以几何作图的准则。
　　高斯还将复数引进了数论，开创了复整数算术理论，复整数在高斯以前只是直观地被引进。
　　高斯是最早怀疑欧几里得几何学是自然界和思想中所固有的那些人之一。
　　当1830年前后匈牙利的波尔约(Janos Bolyai)和俄国的罗巴切夫斯基独立地发表非欧几何学时，高斯宣称他大约在30年前就得到同样的结论。
　　1830年前后，极值(极大和极小)原理在高斯的物理问题和数学研究中开始占有重要地位，例如流体保持静止的条件等问题。
　　在探讨毛细作用时，他提出了一个数学公式能将流体系统中一切粒子的相互作用、引力以及流体粒子和与它接触的固体或流体粒子之间的相互作用都考虑在内。
　天文：
　　1801天文界正在为火星和木星间庞大的间隙烦恼不已，认为火星和木星间应该还有行星未被发现。
　　他在《天体运动理论》(1809)中叙述的方法今天仍在使用，只要稍作修改就能适应现代计算机的要求。高斯在小行星「智神星」方面也获得类似的成功。
　　考虑到其他行星对智神星轨道的摄动，高斯改进了他的计算。这时他的声名远播，荣誉滚滚而来。自那以后，行星、大行星（海王星）接二连三地被发现了。
　　在1812年，他研究了超几何级数，并且把研究结果写成专题论文，呈给哥廷根皇家科学院。
　地理：
　　1820年前后，高斯把注意力转向大地测量——用数学方法测定地球表面的形状和大小。
　　为了增加测量的精确度，他发明了回光仪(一种利用日光以保证比较精确测量的仪器)。
　　他还引进了所谓的高斯误差曲线，并指出概率如何能用变差的钟形曲线(一般称为正态曲线，它是刻画数据统计分布的基础)来表示。
　　1820到1830年间，高斯为了测绘汗诺华公国的地图，开始做测地的工作，他写了关于测地学的书，由于测地上的需要，他发明了日观测仪。
　个人著作：
　　《算术研究》1801年介绍了同余、二次互逆定理
　　《天体运动理论》1809年
　　天体运动的著作
　　《曲面的一般研究》1827年阐述了空间曲面的微积分几何学
　　关于代数基本定理的博士论文1799年证明了每个复系数方程必有复数解
　　《高等大地测量学理论》上1843/44年地理测量
　　《高等大地测量学理论》下1846/47年地理测量
　　《地磁的一般理论》1839年 
　　《地磁概念》1840年 
　　《论与距离平方成反比的引力和斥力的普遍定律》1840年

 

 

 

　　不扯了，go on。

　　就拿luogu 3389【模板】高斯消元法 为例来讲起。

　　题目样例给出：

　　1 3 4 5

　　1 4 7 3

　　9 3 2 2

以此我们可以来列不等式组：

　　X₁ + 3X_{2 + 4X3} = 5.

　　X₁ + 4X₂ + 7X₃ = 3.

　　9X₁ + 3X₂ + 2X₃ = 2.

接下来就是进行消元了，程序当然没人脑辣么灵活，它可不会自动按照最好的方式（加、减）选择。

所以我们需要自己选择一种方式，我习惯用两式相减化简。

步骤：

　　注意这里用a[i][n+1]表示的i行的方程组的值，像上面的5、3、2 .

　　1.定义一个极小值

　　在编程时是存在精度误差的double也只有小数点后14位。

　　判断某个数是否等于0时，小数点后位数少还可以，若是0.0000000000001 呢？是否等于0 呢？答案显然。

　　所以我们需要自定义一个函数。

int sign(double x)            //这是一个极小值
{
    if(fabs(x)<=eps)return 0;    //与0极为接近直接判为0.
    if(x>0)return 1;        //不解释 
    return -1;
}

 

　　2.方程组无解：

　　消元过程中，细心地同学早已发现第i行的方程组里还剩下n-i+1个未知量，另一方面说a[i][i]不等于0，否则无解。

　　所以在解方程时，我们需要特判下。

        if(sign(a[i][i])==0)
        {
            printf("No Solution\n");
            return 0;
        }

 

　　3.两式相减我们需要让方程组上侧式子的系数大于其下方程组，后进行数替换。

        int p=i;
        for(int j=i;j<=n;j++)
        {
            if(fabs(a[j][i])>fabs(a[p][i]))p=j;
        }
        for(int j=1;j<=n+1;j++)swap(a[i][j],a[p][j]);

　　4.然后进行每一步的加减，也就是消元的过程。

        for(int j=i+1;j<=n;j++)
        {
            double ratio=a[j][i]/a[i][i];
            for(int k=1;k<=n+1;k++)a[j][k]=a[j][k]-a[i][k]*ratio;
        }     

 

　　每一步运行情况如下：

　　①

　　②

　　③

　　5.最后一步，当然是求解了，上边过程中你会发现，最后一行只有一个未知数，那么我们可以从上往下一步一步将未知数全解出来。

    for(int i=n;i>=1;i--)
    {
        for(int j=n;j>i;j--)
        {
            if(sign(a[i][j])==0)break;    //找到当前行还不知道准确值得未知数。 
            a[i][n+1]-=ans[j]*a[i][j];    //将所有已知变量的和，从总和中减去。 
        }
        ans[i]=a[i][n+1]/a[i][i];        //用它的总和除以它的系数得变量值，储存答案。 
    }

　　

　　6.输出答案，perfect ！！！

完整代码：

#include <iostream>
#include <cstdio>
#include <cmath>
using namespace std;
int n;
double a[111][111];
double b[111];
double ans[111];
const double eps=1e-8;
int sign(double x)            //这是一个极小值
{
    if(fabs(x)<=eps)return 0;    //与0极为接近直接判为0.
    if(x>0)return 1;        //不解释 
    return -1;
}
int main()
{
//    printf("%lf",eps);
    scanf("%d",&n);
    for(int i=1;i<=n;i++)
    {
        for(int j=1;j<=n+1;j++)scanf("%lf",&a[i][j]);
    }
//    cout<<n<<"\n";
    for(int i=1;i<=n;i++)
    {
        if(sign(a[i][i])==0)
        {
            printf("No Solution\n");
            return 0;
        }
        int p=i;
        for(int j=i;j<=n;j++)
        {
            if(fabs(a[j][i])>fabs(a[p][i]))p=j;
        }
        for(int j=1;j<=n+1;j++)swap(a[i][j],a[p][j]);
//        swap(b[i],b[p]);
        for(int j=i+1;j<=n;j++)
        {
            double ratio=a[j][i]/a[i][i];
            for(int k=1;k<=n+1;k++)a[j][k]=a[j][k]-a[i][k]*ratio;
        }     
//        for(int j=1;j<=n;j++)printf("%lfX1+%lfX2+%lfX3=%lf\n",a[j][1],a[j][2],a[j][3],a[j][4]);cout<<"\n";
    }
    for(int i=n;i>=1;i--)
    {
        for(int j=n;j>i;j--)
        {
            if(sign(a[i][j])==0)break;    //找到当前行还不知道准确值得未知数。 
            a[i][n+1]-=ans[j]*a[i][j];    //将所有已知变量的和，从总和中减去。 
        }
        ans[i]=a[i][n+1]/a[i][i];        //用它的总和除以它的系数得变量值，储存答案。 
    }
    for(int i=1;i<=n;i++)printf("%.2lf\n",ans[i]);
}

（完）

作者：RMY

出处：https://www.cnblogs.com/rmy020718/p/9192236.html