C语言求解线性方程组

2018-02-10 13:07:22来源:cnblogs.com作者:Kukucao人点击

分享

经典问题用高斯约当算法求解线性方程组。这里要求对任意形式的线性方程组都能够妥善处理,不能只适用于方程个数和未知量数目相等的特殊情形。

  先用循环结构将增广矩阵转换为阶梯形矩阵,循环结束时得到阶梯型矩阵非零行行数,同时得到一个链表其中存放有各非零行主元的列标,列标在链表中按从左到右的顺序依次递减。然后根据线性代数中线性方程组的解的情况及判别准则判断方程是否有解,有多少个解。当线性方程组有解时,需要用convert函数将其转换为简化行阶梯型矩阵,然后输出唯一解或一般解

C语言代码如下:

  1 #include <stdio.h>  2 #include <malloc.h>  3 #include <math.h>  4 #define N 5  //增广矩阵列数  5 #define M 3  //增广矩阵行数  6 struct maincol  7 {  8     int col;       //存放各主元下标的结构体类型  9     struct maincol *next; 10 }; 11 typedef struct maincol mc1; 12  13 int test(int s, int t, float a[][N]);  //判断增广矩阵的s行至M行与t列至N列相交形成的子矩阵是否为零矩阵,若是返回0,若不是返回第一个不为零的列的列标 14 void add(mc1 *head, int col, mc1 **tail);   //函数,用于新建一节点,其中保存有主元col列标,然后按递减顺序将其插入maincol类型的链表 15 void convert(float a[][N], int row, mc1 *head);  //函数,用于将阶梯型矩阵转化为简化行阶梯形矩阵 16  17 void main() 18 { 19     float a[M][N];   //增广矩阵 20     char str[N+1]; 21     int i, j; 22     int s, t;    //子矩阵左上角元素行列下标 23     int row, col;  //row用于存放阶梯形矩阵非零行行数 24     float swap; 25     mc1 *head, *tail, *psnew; 26  27     for (i=0; i<M; i++)     //输入并初始化增广矩阵 28     { 29         printf("请输入增广矩阵第%d行/n", i+1); 30         scanf("%s", str); 31         for (j=0; j<N; j++) 32             a[i][j]=str[j]-48; 33     } 34  35     head=(mc1 *) malloc(sizeof(mc1)); 36     head->next=NULL; 37     tail=head; 38     s=t=1;  //子矩阵即为增广矩阵本身,用增广矩阵左上角元素行列标初始化s,t 39  40     while((col=test(s, t, a))!=0)  //子矩阵不为零矩阵 41     { 42         if (s==M)  //增广矩阵已化为阶梯形矩阵 43         { 44             row=s;  //记录非零行行数 45             add(head, col, &tail);  //最后一个非零行主元列标放入maincol类型链表 46             break;     //结束循环 47         } 48         else 49         { 50             j=s-1; 51             for (i=s; i<M; i++) 52             { 53                 if (fabs(a[j][col-1])<fabs(a[i][col-1]))    //列选主元 54                     j=i; 55             } 56  57             if (s-1!=j) 58             { 59                 for (i=col-1; i<N; i++) 60                 { 61                     swap=a[j][i]; 62                     a[j][i]=a[s-1][i];    //列选主元 63                     a[s-1][i]=swap; 64                 } 65             } 66  67             if (col==N)     //增广矩阵已经化为阶梯形矩阵 68             { 69                 row=s;         //记录非零行行数 70                 add(head, col, &tail);      //最后一个非零行主元列标放入maincol类型链表 71                 break;    //结束循环 72             } 73  74             for (i=s; i<M; i++) 75                 a[i][col-1]=-(a[i][col-1]/a[s-1][col-1]); 76  77             for (i=col; i<N; i++)                        //消元 78             { 79                 for (j=s; j<M; j++) 80                     a[j][i]=a[j][col-1]*a[s-1][i]+a[j][i]; 81             } 82  83             add(head, col, &tail);      //将消元后得到的主元列标col放入maincol类型链表 84             s++; 85             t=col+1;          //更新s,t,使s,t成为消元后得到的新的子矩阵的左上角元素行列标,为test函数操作新的子矩阵作准备 86             continue;      //开始新一轮循环 87         } 88     } 89     if (col==0)    //从循环控制条件退出循环 90         row=s-1;   //此时增广矩阵已成为阶梯形矩阵,非零行函数就是s-1 91  92     if (row==0)          //利用线性方程组解的判别准则判断是否有解,有多少个解 93     { 94         printf("线性方程组有无穷多组解/n");  //增广矩阵为零矩阵,无穷多组解 95         printf("一般解:/n"); 96         for (i=1; i<N; i++) 97             printf("x%d=t%d/n", i, i);   //输出解 98     } 99     else100     {101         psnew=head->next;102         if (psnew->col==N)       //阶梯形矩阵最后一主元在最后一列,无解103             printf("线性方程组无解/n");                104         else105         {106             convert(a, row, head);   //用convert函数将阶梯形矩阵进一步化为简化行阶梯形矩阵107             if (row==N-1)    //非零行行数等于未知量个数,有唯一解108             {109                 printf("线性方程组有唯一解:/n");110                 for (i=1; i<=row; i++)                   //输出唯一解111                     printf("x%d=%f/n", i, a[i-1][N-1]);112             }113             else    //非零行行数小于未知量个数,有无穷多组解114             {115                 int *m;116                 m=(int *) malloc((N-1-row)*sizeof(int));117 118                 i=N-1-row;119                 for (j=N-1; j>=1; j--)120                 {121                     if (j!=psnew->col)122                     {123                         m[--i]=j;        //从所有未知量标号中筛选出自由未知量标号           124                         if (i==0)125                             break;126                     }127                     else128                     {129                         if (psnew->next!=NULL)130                             psnew=psnew->next;131                     }132                 }133                 printf("线性方程组有无穷多组解/n");134                 printf("一般解:/n");135                 i=row;136                 for (psnew=head->next; psnew!=NULL; psnew=psnew->next)137                 {138                     printf("x%d=%f", psnew->col, a[i-1][N-1]);     //输出一般解139                     for (j=0; j<N-1-row; j++)140                     {141                         if(m[j]<psnew->col)142                         {143                             printf("-%dx%d", 0, m[j]);144                         }145                         else146                         {147                             printf("-%fx%d", a[i-1][m[j]-1], m[j]);148                         }149                     }150                     printf("/n");151                     i--;152                 }153             }154         }155     }156 }157 158 int test(int s, int t, float a[][N])        //判断增广矩阵的s行至M行与t列至N列相交形成的子矩阵是否为零矩阵,若是返回0,若不是返回第一个不为零的列的列标159 {160     int i, j;161 162     for (j=t-1; j<N; j++)163     {164         for (i=s-1; i<M; i++)165         {166             if (a[i][j]!=0)167                 return (j+1);168         }169     }170     return (0);171 }172 173 void add(mc1 *head, int col, mc1 **tail)    //函数,用于新建一节点,其中保存有主元col列标,然后按递减顺序将其插入maincol类型的链表174 {175     mc1* psnew;176 177     psnew=(mc1 *) malloc(sizeof(mc1));178     psnew->col=col;179 180     if(head->next==NULL)181     {182         psnew->next=NULL;183         head->next=psnew;184         *tail=psnew;185     }186     else187     {188         psnew->next=head->next;189         head->next=psnew;190     }191 }192 193 void convert(float a[][N], int row, mc1 *head)    //函数,用于将阶梯型矩阵转化为简化行阶梯形矩阵194 {195     mc1 *psnew, *pq;196     int i, j, k, m;197 198     psnew=head->next;199     for (i=row-1; i>=0; i--)200     {201         if (a[i][psnew->col-1]!=1)     //各非零行主元化为1202         {203            for (j=psnew->col; j<N; j++)204                a[i][j]=a[i][j]/a[i][psnew->col-1];205         }206 207         psnew=psnew->next;208     }209 210     psnew=head->next;           //向上消元把除第一个主元外其余主元所在列中在主元上方的部分变为零211     for (i=row-1; i>=1; i--)212     {213         m=N-psnew->col-(row-i);  //获取未知量标号1,2,--,N-1中位于i+1非零行主元列标号右边的自由未知量标号个数214         for (j=i-1; j>=0; j--)215         {216             pq=head->next;     //pq指向存放最后一个非零行主元列标号的节点217             for (k=N; k>psnew->col; k--)218             {219                 if (k!=pq->col)220                 {221                    a[j][k-1]=-(a[i][k-1]*a[j][psnew->col-1])+a[j][k-1];  //从右向左作初等行变换直至i+1行主元所在列右边的列位置,期间跳过i+2----row行主元所在的列222                    m--;223                    if (m==0)224                        break;225                 }226                 else227                 {228                     if (pq->next!=psnew)229                         pq=pq->next;230                 }231             }232         }233         psnew=psnew->next;  //递进至上一行主元,为新一轮向上消元作准备234     }235 }

最新文章

123

最新摄影

微信扫一扫

第七城市微信公众平台