高斯消元法求线性方程组

高斯消元法

作用
可以快速求解n元线性方程组：

{\begin{cases} a_{11} x_{1} + a_{12} x_{2} + a_{13} x_{3} + \dots + a_{1 n} x_{n} = b_{1} \\ a_{21} x_{1} + a_{22} x_{2} + a_{23} x_{3} + \dots + a_{2 n} x_{n} = b_{2} \\ \dots \\ a_{n 1} x_{1} + a_{n 2} x_{2} + a_{n 3} x_{3} + \dots + a_{n n} x_{n} = b_{n} \end{cases}

思路
利用线性代数知识，存下方程的增广矩阵

(\begin{array}{ccc} a_{11} & a_{12} & a_{13} & \dots & a_{1 n} & b_{1} \\ a_{21} & a_{22} & a_{23} & \dots & a_{2 n} & b_{2} \\ ⋮ & ⋮ & ⋮ & ⋮ & ⋮ & ⋮ \\ a_{n 1} & a_{n 2} & a_{n 3} & \dots & a_{n n} & b_{n} \end{array})

通过初等行变换将其化为行阶梯型矩阵，即可求解。
行变换步骤：
（1）对每一列进行遍历找到绝对值最大的元素所在的行。^[1]
（2）将这一行的所有元素与当前的第一行所有元素进行交换。
（3）将当前第一行的第一个元素化为1。
（4）从第二行开始把当前第一列的元素全部消为0。

代码实现

#include<iostream>
#include<algorithm>
#include<cmath>
const int N=110;
using namespace std;
const double eps=1e-6;//代表误差
int n;
double a[N][N];
int gauss()
{
    int c,r;//分别表示当前最小列，行
    for(c=0,r=0;c<n;c++)
    {
        int t=r;
        for(int i=r;i<n;i++ )
        {
            if(fabs(a[i][c])>fabs(a[t][c]))//注意浮点数是有误差的，不能直接写等于0。
            {
                t=i;
            }
        }
        if(fabs(a[t][c])<eps)
        {
            continue;
        }
        for(int j=c;j<=n;j++)
        {
            swap(a[t][j],a[r][j]);
        }
        for(int j=n;j>=c;j--)
        {
            a[r][j]/=a[r][c];
        }
        for(int j=r+1;j<n;j++)
        {
            if(fabs(a[j][c])>eps)
            {
                for(int k=n;k>=c;k--)
             {
                a[j][k]-=a[j][c]*a[r][k];
             }
            }
           
        }
        r++;
    }//其实每一次没有被continue的循环就是删去系数矩阵的一条长和一条宽，变为一个较小的正方形。
    if(r<n)//n-r等于几那么最后这几行就全为0
    {
        for(int i=r;i<n;i++)
        {
            if(fabs(a[i][n])>eps)
            {
                return 2;
            }
        }
        return 1;
    }
    for(int i=n-2;i>=0;i--)//从下往上消去每个元素上面的元素
    {
        for(int j=i+1;j<n;j++)
        {
            a[i][n]-=a[j][n]*a[i][j];
        }
    }
     return 0;
}
int main()
{
    cin>>n;
    for(int i=0;i<n;i++)
        {
            for(int j=0;j<=n;j++)
            {
                scanf("%lf",&a[i][j]);
            }
        }
    int res=gauss();
    if(res==0)
    {
        for(int i=0;i<n;i++)
        {
            printf("%.2lf\n",a[i][n]);
        }
    }
    else if(res==2)
    {
        printf("No solution");
    }
    else if(res==1)
    {
        printf("Infinite group solutions");
    }
    return 0;
}