P3389 【模板】高斯消元法

· · 题解

高斯消元其实脱去华丽的外表，便是初中的普通的消元法。就是解多元一次方程的一个算法罢了。例如：

这便是一个极其普通的二元一次方程。以一个普通中学生的身份我们怎么做这题呢？

显而易见，我们将先①3-② 2得9y-10y=-3，解得y=3,将y=3代入①式中便可得x=2。

怎么样？是不是很简单？其实这就是高斯消元得核心思想。高斯消元不过是将此推广到多元。

例如，如果我们要解以下方程：

我们可以将其系数写成矩阵的形式如：

| 2   1   -3 |

| 1   -1  2 |

| 5   2   1 |

按照我们前面的思路，我们将第二行每项系数乘上2/1，即一式和二式第一项系数之比，然后将其相减，存在第二行，那么第二行便成为了| 0 3 -7 |。每行都这么操作，便能将除了第一行以外的第一项系数清0。得：

| 2   1   -3 |

| 0   3   -7 |

| 0   1/5  -17/5|

这样就能消去第一项，以此类推，便能将（n-1）个系数消去得到一个上三角的矩阵。

最后便能解得：x=2,y=5,z=-4。

最后的最后，贴上我丑陋的代码（溜）：

#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;

int n;
double ans[110];  //ans数组存放答案 
double fc[110][110];

void Gauss_Work1()  //消元 
{
    int i,j,k,zx;
    double p;
    for(i=1;i<=n;i++)
    {
        for(zx=i;zx<=n;zx++)
            if(fc[zx][i]!=(double)0)  //先去寻找该项系数不为0的方程式，否则没意义 
                break;
        if(zx==n+1)
            <%printf("No Solution");exit(0);%>  //如果没有找到该项系数不为0的方程式，则说明此项未知数取任何数都没有关系，称为自由元 
        for(j=i;j<=n+1;j++)
            swap(fc[i][j],fc[zx][j]);  //将找到的方程式与第i行交换，便保证了第i行的方程式此项非0，能用此式消去其它方程式的该项 
        for(j=i+1;j<=n;j++)
        {
            if(fc[j][i]==0)  //如果此项为0，则没有必要进行消元 
                continue;
            p=fc[i][i]/fc[j][i];   
            for(k=i;k<=n+1;k++)
                fc[j][k]=fc[i][k]-fc[j][k]*p;  //将两式相减，赋值于第j行，成功消去此元 
        }
    }
}

void Gauss_Work2()  //计算各项未知数的值 
{
    int i,j;
    for(i=n;i>=1;i--)
    {
        ans[i]=fc[i][n+1];  //先将ans数组赋值为消元后方程的值 
        for(j=i+1;j<=n;j++)
            ans[i]-=fc[i][j]*ans[j];  //将方程式减去各项已知未知数和系数的乘积 
        ans[i]/=fc[i][i];  //除于你要求的未知数所对应的系数，即解类似kx=y(k,y已知)的方程式 
    }
    //算法结束(撒花) 
}

int main()
{
    int i,j;
    scanf("%d",&n);  //输入n
    for(i=1;i<=n;i++)
        for(j=1;j<=n+1;j++)
            scanf("%lf",&fc[i][j]);  //输入方程式
    Gauss_Work1();  //消元 
    Gauss_Work2();  //计算各项未知数的值 
    for(i=1;i<=n;i++)
        printf("%.2lf\n",ans[i]);  //输出答案 
    return 0;
}

对了，记得变量是双精度浮点型的(废话)。