高斯消元法，最后两位小数不正确，估计是算法有问题，但具体找不出来

Question

问题遇到的现象和发生背景

windows10操作系统，visual studio2019编译器

问题相关代码，请勿粘贴截图

已经在函数内定义了方程组的各个参数，要求高斯消元法解方程组
原代码：

#include<stdio.h>
#include<math.h>
int main()
{
    long double a[4][4] = { 1.1161,0.1254,0.1397,0.1490,0.1582,1.1675,0.1768,0.1871,0.2368,0.2471,0.2568,1.2671,0.1968,0.2071,1.2168,0.2271 };
    long double b[4] = { 1.5471,1.6471,1.8471,1.7471 };
    long double change[4], exchange;
    long double ans[4];
    long double right,temp;
    long double max = 0;
    int i, j, k;
    void print();
    print();
    //列选主元
    for (k = 0;k <= 3;k++)
    {
        max = 0;
        for (i = k;i <= 3;i++)
            max = max < a[i][k] ? a[i][k] : max;
        for(i=k;i<=3;i++)
            if (a[i][k]==max)
            {
                for (j = 0;j <= 3;j++)
                {
                    change[j] = a[k][j];
                    a[k][j] = a[i][j];
                    a[i][j] = change[j];
                    exchange = b[k];
                    b[k] = b[i];
                    b[i] = exchange;
                }
            }
    //转换成阶梯型矩阵
        for (i = k + 1;i <= 3;i++)
        {
            temp = a[i][k] / a[k][k];
            for (j = k;j <= 3;j++)
            {
                a[i][j] = a[i][j] - a[k][j] * temp;
            }
            b[i] = b[i] - b[k] * temp;
            
        }
    }
    //求回代求解
    for (k = 3;k >= 0;k--)
    {
        right = b[k];
        for (i = k + 1;i <= 3;i++)
            right =right - a[k][i] * ans[i];
        ans[k] = right / a[k][k];
    }
    printf("\nThe ans =\n");
    for (k = 0;k <= 3;k++)
        printf("\tx(%d)=%.4lf\n", k + 1, ans[k]);
    return 0;
}
//函数，用来输出系数矩阵A和向量B
void print()                 
{
    double a[4][4] = { 1.1161,0.1254,0.1397,0.1490,0.1582,1.1675,0.1768,0.1871,0.2368,0.2471,0.2568,1.2671,0.1968,0.2071,1.2168,0.2271 };
    double b[4] = { 1.5471,1.6471,1.8471,1.7471 };
    int n,j;
    printf("MAT A = \n");
    for (n = 0;n <= 3;n++)
    {
        for (j = 0;j <= 3;j++)
        {
            printf("\t%.4lf", a[n][j]);
            if (j == 3)printf("\n");
        }
    }
    printf("MAT B = \n");
    for (n = 0;n <= 3;n++)
        printf("\t%.4lf\n", b[n]);
}                 

运行结果及报错内容

我想要达到的结果

希望最后两个x3为0.9351，x4为0.881左右

Answer 1

不需要列选主元，按照正常高斯消元法步骤，将增广矩阵转为阶梯型矩阵后，回代求解即可，把列选主元的步骤删除就得到你想要的答案了

#include<stdio.h>
#include<math.h>
int main()
{
    long double a[4][4] = { 1.1161,0.1254,0.1397,0.1490,0.1582,1.1675,0.1768,0.1871,0.2368,0.2471,0.2568,1.2671,0.1968,0.2071,1.2168,0.2271 };
    long double b[4] = { 1.5471,1.6471,1.8471,1.7471 };
    long double change[4], exchange;
    long double ans[4];
    long double right,temp;
    long double max = 0;
    int i, j, k;
    void print();
    print();
    for (k = 0;k <= 3;k++)
    {
    //转换成阶梯型矩阵
        for (i = k + 1;i <= 3;i++)
        {
            temp = a[i][k] / a[k][k];
            for (j = k;j <= 3;j++)
            {
                a[i][j] = a[i][j] - a[k][j] * temp;
            }
            b[i] = b[i] - b[k] * temp;
        }
    }
    //求回代求解
    for (k = 3;k >= 0;k--)
    {
        right = b[k];
        for (i = k + 1;i <= 3;i++)
            right =right - a[k][i] * ans[i];
        ans[k] = right / a[k][k];
    }
    printf("\nThe ans =\n");
    for (k = 0;k <= 3;k++)
        printf("\tx(%d)=%.4Lf\n", k + 1, ans[k]);
    return 0;
}
//函数，用来输出系数矩阵A和向量B
void print()                 
{
    double a[4][4] = { 1.1161,0.1254,0.1397,0.1490,0.1582,1.1675,0.1768,0.1871,0.2368,0.2471,0.2568,1.2671,0.1968,0.2071,1.2168,0.2271 };
    double b[4] = { 1.5471,1.6471,1.8471,1.7471 };
    int n,j;
    printf("MAT A = \n");
    for (n = 0;n <= 3;n++)
    {
        for (j = 0;j <= 3;j++)
        {
            printf("\t%.4Lf", a[n][j]);
            if (j == 3)printf("\n");
        }
    }
    printf("MAT B = \n");
    for (n = 0;n <= 3;n++)
        printf("\t%.4Lf\n", b[n]);
}