| 网站首页 | 业界新闻 | 小组 | 威客 | 人才 | 下载频道 | 博客 | 代码贴 | 在线编程 | 编程论坛
欢迎加入我们,一同切磋技术
用户名:   
 
密 码:  
共有 638 人关注过本帖
标题:这个高斯列主元消去法,哪里错了?
只看楼主 加入收藏
dhz662820909
Rank: 1
等 级:新手上路
帖 子:12
专家分:0
注 册:2013-10-28
结帖率:40%
收藏
已结贴  问题点数:19 回复次数:4 
这个高斯列主元消去法,哪里错了?
#include<stdio.h>
#include<stdlib.h>
#include<math.h>
#define N 4
#define Error 0.00001
void Gauss(float p[N][N+1],float *p1);
void Solve(float p[N][N+1],float *p1,float *p2,float *p3,int rank,int a_rank);
void Print(float p[N][N+1]);
void Solve_2(float p[N][N+1],float *p1,int rank,float *p3);
int main(void)
{
    int i,j;
    float p[N][N+1],*p1=NULL;     //用于调用高斯列主元消去法后,返回并打印方程唯一解
    p1=(float*)calloc(N,sizeof(float));
    if(p1==NULL)
    {
                printf("No enough memory(main_p)!\n");
                exit(1);
    }
    for(i=0;i<N;i++)     //输入增广矩阵
    {
                    printf("Row %d,enter %d elements:\n",i+1,N+1);
                    for(j=0;j<N+1;j++)
                    {
                                    scanf("%f",&p[i][j]);
                    }
    }
    Gauss(p,p1);   //调用高斯列主元消去法
    printf("\n");
    free(p1);
    system("pause");
    return(0);
}
void Gauss(float p[N][N+1],float *p1)       //p1用于保存解
{
     int i,j,k,index,revise=0,s,t,rank=0,a_rank=0;   //index是某一列主元的下标,revise起到修正主元下标的作用
     float *p2=NULL,*p3=NULL,Coeff;    //p2用于换行以及打印基础解系,p3用于保存主元及其所在的列
     p2=(float*)calloc(N,sizeof(float));
     p3=(float*)calloc(N,sizeof(float));
     if(p2==NULL)
     {
                 printf("No enough memory(Gauss_p2)!\n");
                 exit(1);
     }
     if(p3==NULL)
     {
                 printf("No enough memory(Gauss_p3)!\n");
                 exit(1);
     }
     for(j=0;j<N;j++)    //表示第j列选主元
     {
                     index=j-revise;
                     for(i=j+1-revise;i<N;i++)   //找到第j列主元的下标位置
                     {
                                    
                                     if(fabs(p[index][j])+Error<fabs(p[i][j]))
                                     index=i;
                     }
                     for(k=j;k<N+1;k++)   //第j-revise行与主元所在的行交换
                     {
                                     p2[k]=p[j-revise][k];
                                     p[j-revise][k]=p[index][k];
                                     p[index][k]=p2[k];
                     }
                     if(fabs(p[j-revise][j])<Error)      //p[j-revise][j]是主元
                     {
                                                   revise=revise+1;
                                                   continue;
                     }
                     else
                     {
                         rank=rank+1;           //计算系数矩阵的秩
                         p3[j]=p[j-revise][j];   //把主元保存在p3数组中
                     }
                     for(s=j+1-revise;s<N;s++)    //消元
                     {                 
                                       Coeff=p[s][j]/p[j-revise][j];
                                       for(t=j;t<N+1;t++)
                                       p[s][t]=p[s][t]-Coeff*p[j-revise][t];
                     }
     }
     /*printf("\nPrimary:\n");   //打印p3数组,显示主元和主元所在的列
     for(i=0;i<N;i++)
     printf("%f ",p3[i]);
     printf("\n");*/
     for(i=N;i>=0;i--)   //计算增广矩阵的秩
     {
                      if(fabs(p[i][N])>Error)
                      {
                                             a_rank=i+1;
                                             break;
                      }
     }
     Solve(p,p1,p2,p3,rank,a_rank);
     free(p2);
     free(p3);
}
void Solve(float p[N][N+1],float *p1,float *p2,float *p3,int rank,int a_rank)
{
     int i,j,k;
     float term,sum;
     if(rank==N)   //系数矩阵满秩 ,方程有唯一解
     {
                   printf("\nAfter Gauss:\n");   //打印经过高斯列主元消去法后的增广矩阵
                   Print(p);
                   for(i=N-1;i>=0;i--)
                   {
                                      sum=0;
                                      for(j=N-1;j>i;j--)
                                      {
                                                          term=p[i][j]*p1[j];
                                                          sum=sum+term;
                                      }
                                      p1[i]=(p[i][N]-sum)/p3[i];
                   }
                   printf("\nSolution:");     //打印方程的唯一解
                   for(i=0;i<N;i++)
                   printf("%f ",p1[i]);
                   printf("\n");
                  
     }
     else if(rank>=a_rank && rank<N)  //系数矩阵不满秩,系数矩阵的秩等于增广矩阵的秩,方程有无穷多解
     {
          printf("\nAfter Gauss:\n");   //打印经过高斯列主元消去法后的增广矩阵
          Print(p);
          Solve_2(p,p1,rank-1,p3);
          for(i=N-1;i>=0;i--)
          {
                             sum=0;
                             for(j=N-1;j>i;j--)
                             {
                                               term=p[i][j]*p1[j];
                                               sum=sum+term;
                             }
                             if(fabs(p3[i])>Error)
                             p1[i]=(p[i][N]-sum)/p3[i];
          }
          printf("\nParticular solution:\n");
          for(i=0;i<N;i++)
          printf("%f ",p1[i]);
          for(i=0;i<N;i++)
          p[i][N]=0;
          for(k=N-1;k>=0;k--)
          {
                             if(fabs(p3[k])<Error)
                             {
                                                  p1[k]=1;
                                                  for(i=N-1;i>=0;i--)
                                                  {
                                                                     sum=0;
                                                                     for(j=N-1;j>i;j--)
                                                                     {
                                                                                       term=p[i][j]*p1[j];
                                                                                       sum=sum+term;
                                                                     }
                                                                     if(fabs(p3[i])>Error)
                                                                     p1[i]=(p[i][N]-sum)/p3[i];
                                                  }
                             printf("\nGeneral solution:\n");
                             for(i=0;i<N;i++)
                             printf("%f ",p1[i]);
                             p1[k]=0;
                             }
          }
     }
     else if(rank<a_rank)    //系数矩阵的秩不等于增广矩阵的秩,方程无解
     {
          printf("\nAfter Gauss:\n");   //打印经过高斯列主元消去法后的增广矩阵
          Print(p);
          printf("\nNo solution!\n");
     }
}
void Print(float p[N][N+1])
{
     int i,j;
     for(i=0;i<N;i++)
     {
                     for(j=0;j<N+1;j++)
                     printf("%f ",p[i][j]);
                     printf("\n");
     }
}
void Solve_2(float p[N][N+1],float *p1,int rank,float *p3)     //p1和p3就是Solve中的p1和p3
{
     int i,j,k,s,t;
     float *p2=NULL,*p4=NULL;
     p2=(float*)calloc(N,sizeof(float));       //用来换行
     p4=(float*)calloc(N,sizeof(float));       //用来保存p3
     if(p2==NULL)
     {
                 printf("No enough memory(Solve_p2)!\n");
                 exit(1);
     }
     if(p4==NULL)
     {
                 printf("No enough memory(Solve_p4)!\n");
                 exit(1);
     }
     for(i=0;i<N;i++)
     p4[i]=p3[i];
     for(i=rank;i>=0;i--)
     {
                           for(j=N-1;j>=0;j--)
                           {
                                              if(fabs(p4[j])>Error)
                                              {
                                                                   for(k=j;k<N+1;k++)
                                                                   {
                                                                                     p2[k]=p[j][k];
                                                                                     p[j][k]=p[i][k];
                                                                                     p[i][k]=p2[k];
                                                                   }
                                                                   p4[j]=0;
                                                                   break;
                                              }
                           }
     }
     free(p2);
     free(p4);
}
这个是用Dev C++写的,同样一组数据,有时候可以运行,有时候不可以。。。
搜索更多相关主题的帖子: 高斯 include 
2013-12-26 12:16
so_love
Rank: 13Rank: 13Rank: 13Rank: 13
等 级:蒙面侠
威 望:7
帖 子:812
专家分:4151
注 册:2013-11-25
收藏
得分:7 
额。代码好长。排版有点乱。 而且不太明白什么叫有时好使有时不好使?还是同一组数据。你这有具体模块。找相应的模块呗

一花一世界、一叶一追寻、片片花叶落、情系何人身。
2013-12-26 12:40
xiaozi2013
Rank: 13Rank: 13Rank: 13Rank: 13
来 自:四川成都
等 级:蒙面侠
威 望:6
帖 子:830
专家分:4331
注 册:2013-12-13
收藏
得分:7 
-------最近写代码,特别在意排版,及功能模块的区分-----------
个人觉得:

-------还不错哦------------

--整天敲那些破代码能找着媳妇儿吗。。。---
2013-12-26 13:20
so_love
Rank: 13Rank: 13Rank: 13Rank: 13
等 级:蒙面侠
威 望:7
帖 子:812
专家分:4151
注 册:2013-11-25
收藏
得分:0 
这个是有排版。但是距离太大了。看着不舒服

一花一世界、一叶一追寻、片片花叶落、情系何人身。
2013-12-26 17:08
beyondyf
Rank: 19Rank: 19Rank: 19Rank: 19Rank: 19Rank: 19
等 级:贵宾
威 望:103
帖 子:3282
专家分:12654
注 册:2008-1-21
收藏
得分:7 
唉,太让我失望了。什么乱七八糟的东西,我耐着性子看完了Gauss和Solve两个函数,各种废话和错误。

重剑无锋,大巧不工
2013-12-26 18:39
快速回复:这个高斯列主元消去法,哪里错了?
数据加载中...
 
   



关于我们 | 广告合作 | 编程中国 | 清除Cookies | TOP | 手机版

编程中国 版权所有,并保留所有权利。
Powered by Discuz, Processed in 0.040158 second(s), 7 queries.
Copyright©2004-2024, BCCN.NET, All Rights Reserved