C语言实现列主元消去法

首先说明一下列主元消去法的内容,以下是360搜索到的,应该可以看懂

 

方法说明(以4阶为例):

第1步消元--在增广矩阵(A,b)第一列中找到绝对值最大的元素,将其所在

行与第一行交换,再对(A,b)做初等行变换使原方程组转化为如下形式:

第2步消元--在增广矩阵(A,b)中的第二列中(从第二行开始)找到绝对值

最大的元素,将其所在行与第二行交换,再对(A,b)做初等行变换使原方程组转化为:

第3步消元--在增广矩阵(A,b)中的第三列中(从第三行开始)找到绝对值

最大的元素,将其所在行与第三行交换,再对(A,b)做初等行变换使原方程组转化为:

按x4 → x3→ x2→ x1 的顺序回代求解出方程组的解。

 

然后C语言实现   参考了一位大神的博客 :   https://blog.csdn.net/u014070279/article/details/48317905

 

//获取主元所在行
int getMaxinumLaber( double a[N][M], int k)
{
    int i;
    int laber=k;
    double maxinum=0;

    //列固定,行变,找最大值
    for(i=k; i<N; i++)
        if(maxinum<fabs(a[i][k]))
        {
            maxinum = fabs(a[i][k]);
            laber = i;
        }
    return laber;
}

//列主元高斯消去法
void gaussian_elimination( double a[N][M],double *x)
{
   int i, j, k;
    int laber;
    double temp;
   // double **a = a_mai;
    //double *x=x_mai;

    //高斯消元法
    for(k=0; k<N; k++)
    {
        laber = getMaxinumLaber(a, k);


        //主元交换,就是laber行和第k行交换
        if(laber!= k)
        {
          for(i=0; i<M; i++)
          {
              temp = a[k][i];
              a[k][i] = a[laber][i];
              a[laber][i] = temp;
          }
        }


        //消元
        for(i=k+1; i<N; i++)
        {

            if(a[k][k]==0)

            break;


            temp = a[i][k]/a[k][k];


            for(j=k; j<M; j++)
            a[i][j] = a[k][j]*temp-a[i][j];
        }
    }


    //回代求解x
    double sum;
    for(i=N-1; i>=0; i--)
    {
        sum = 0;
        for(j=i+1; j<N; j++)
        {
          sum += a[i][j]*x[j];
        }
        x[i] = (a[i][M-1]-sum)/a[i][i];
    }

    printf("\n");
    printf("变换增广系数矩阵\n");
    for(i=0; i<N; i++)
    {
        for(j=0; j<M; j++)
        printf("%f   ", a[i][j]);
        printf("\n");
    }


    printf("\n");
    for(i=0; i<N; i++)
    printf("%f   ", x[i]);
    printf("\n");

}

然后讲一下如何用:

a[N][M]其实就是我们线性代数中的增广系数矩阵,x[N]就是待求解的未知数向量;

应用是结合最小二乘法进行曲线拟合,首先通过最小二乘法求得增广系数矩阵,然后就是求拟合函数y =ax+b中的系数,b就是x[0],a就是x[1]   .

如果不了解最小二乘法的话就直接找个方程试试。

int main(void)
{
    double x_know[5] ={1,2,3,4,5},y_know[5] ={7,11,17,27,40},x[N]={0};
    int n =5,i=0,j=0;
    double A[N][M]={0};

    for(i=0;i<n;i++)
    {

      A[0][0]+=1;
      A[0][1]+=x_know[i];
      A[0][2]+=y_know[i];
      A[1][0]=A[0][1];
      A[1][1]+=x_know[i]*x_know[i];
      A[1][2]+=x_know[i]*y_know[i];
    }

    printf("\n");
    printf("原增广系数矩阵\n");
    for(i=0; i<N; i++)
    {
        for(j=0; j<M; j++)
        printf("%f   ", A[i][j]);
        printf("\n");
    }
    gaussian_elimination(A,x);

    return 0;
}

  • 4
    点赞
  • 20
    收藏
    觉得还不错? 一键收藏
  • 打赏
    打赏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

GallagherZ

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付:¥1
获取中
扫码支付

您的余额不足,请更换扫码支付或充值

打赏作者

实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值