高斯消元法（高斯·约当消元法）（浮点）

高斯消元法，用于解多元一次方程（几乎类似模拟手动解方程）。

思路：
通过等式的乘除，把方程1的$x_1$系数$a_{11}$分别化为方程2~方程n的$x_1$系数，然后将方程2~方程n减去得到的新方程，从而消掉方程2~方程n中的$x_1$。接着用方程2的$x_2$继续把方程3~方程n中的$x_2$消掉……
大概系数就成了这个样子↓

举个例子：

一个问题：遇到这种情况

$$\begin{matrix} 0.0000000000001×1&+&x_2&=&1\\ x_1&+&x_2&=&2\\ \end{matrix}$$
就会出现精度问题，处理成
$$\begin{matrix} x_1&+&1000000000000x_2&=&1000000000000\\ 0&+&999999999999x_2&=&2\\ \end{matrix}$$
而当位数一多，我们的计算机就有可能存不了那么多9，从而忽略掉一堆，原本的误差就扩散得越来越大，这时候我们需要选主消元（比如把刚才两个式子交换一下），选择当前这一项系数绝对值最大的方程作为主元，来消掉其它方程。

高斯·约当消元：
高斯消元是当执行到第i个方程的$x_k$时，消掉i以后的$x_k$。而约当就同时消掉i以前的$x_k$，使系数变成一条对角线

解得情况：高斯消元完成后，下面若有方程系数全部为0
$0×1+0×2+0×3+...+0×n=0$
则说明有多解
如果为
$0×1+0×2+0×3+...+0×n=非0$
则无解

代码：

#include<cstdio>
#include<cmath>
#include<algorithm>
using std::swap;
#define MAXN 405
#define EPS 1e-8
typedef double Matrix[MAXN][MAXN];//系数矩阵
int n,m;
int Rank;//有效方程的行数，Rank之后的方程x系数为0
double X[MAXN];//解
bool Free[MAXN];//是否为自由变量
Matrix A;//系数矩阵
//浮点数比较
int fcmp(double a,double b)
{
    if((a-b)>EPS)
        return 1;
    else if((a-b)>-EPS)
        return 0;
    return -1;
}
//高斯·约当消元
void Gauss()
{
    int r,c,mxr;
    for(r=1,c=1;r<=n&&c<m;r++,c++)
    {
        //寻找绝对值最大（选主）
        mxr=r;
        for(int i=r+1;i<=n;i++)
            if(fcmp(fabs(A[i][c]),fabs(A[mxr][c]))>0)
                mxr=i;
        //第c项在第r个方程之后系数都为0
        if(fcmp(A[mxr][c],0.0)==0)
        {r--;continue;}//执行下一项
        //选好主，交换方程
        if(mxr!=r)swap(A[mxr],A[r]);
        //消元
        for(int i=1;i<=n;i++)
            if(i!=r&&fcmp(A[i][c],0.0)!=0)
                for(int j=m;j>=c;j--)
                    A[i][j]-=A[r][j]/A[r][c]*A[i][c];
    }
    Rank=r-1;
}
//判断是否有解
bool check()
{
    //判断是否无解
    for(int i=Rank+1;i<=n;i++)
        if(fcmp(A[i][m],0.0)!=0)//在x系数为0时，结果不为0
            return 0;
    //初始所有都是未知变量
    for(int i=1;i<m;i++)
        Free[i]=1;
    //计算解
    for(int i=Rank,cnt=0,pos;i>0;i--)
    {
        cnt=0;//记数当前方程未知变量个数
        for(int j=1;j<m;j++)
            if(Free[j]&&fcmp(A[i][j],0.0)!=0)
                cnt++,pos=j;
        if(cnt==1)//一个未知变量，可计算
        {
            Free[pos]=0;
            X[pos]=A[i][m]/A[i][pos];
        }
    }
    return 1;
}
int main()
{
    freopen("Gauss_data.in","r",stdin);
    scanf("%d%d",&n,&m);
    m++;
    for(int i=1;i<=n;i++)
        for(int j=1;j<=m;j++)
            scanf("%lf",&A[i][j]);
    Gauss();
    if(!check())
        printf("No solution\n");
    else
    {
        if(Rank<m-1)
        {
            printf("Multiple solution, free_num: %d\n",m-1-Rank);
            for(int i=1;i<m;i++)
                if(Free[i])
                    printf("X[%d] not determined\n",i);
                else
                    printf("X[%d] = %.4lf\n",i,X[i]);
        }
        else
            for(int i=1;i<m;i++)
                printf("X[%d] = %.4lf\n",i,X[i]);
    }
    return 0;
}

还可以继续看高斯·约当消元法（整数）