C語言求解線性方程組

C語言求解線性方程組

經典問題用高斯約當演算法求解線性方程組。這裡要求對任意形式的線性方程組都能夠妥善處理，不能只適用於方程個數和未知量數目相等的特殊情形。

先用迴圈結構將增廣矩陣轉換為階梯形矩陣，迴圈結束時得到階梯型矩陣非零行行數，同時得到一個鏈表其中存放有各非零行主元的列標，列標在鏈表中按從左至右的順序依次遞減。然後根據線性代數中線性方程組的解的情況及判別準則判斷方程是否有解，有多少個解。當線性方程組有解時，需要用convert函數將其轉換為簡化行階梯型矩陣，然後輸出唯一解或一般解

C語言代碼如下：

#include <stdio.h>
#include <malloc.h>
#include <math.h>
#define N 5  //增廣矩陣列數
#define M 3  //增廣矩陣行數
struct maincol
{
    int col;      //存放各主元下標的結構體類型
    struct maincol *next;
};
typedef struct maincol mc1;

int test(int s, int t, float a[][N]);  //判斷增廣矩陣的s行至M行與t列至N列相交形成的子矩陣是否為零矩陣，若是返回0，若不是返回第一個不為零的列的列標
void add(mc1 *head, int col, mc1 **tail);  //函數，用於建立一節點，其中儲存有主元col列標，然後按遞減順序將其插入maincol類型的鏈表
void convert(float a[][N], int row, mc1 *head);  //函數，用於將階梯型矩陣轉化為簡化行階梯形矩陣

void main()
{
    float a[M][N];  //增廣矩陣
    char str[N+1];
    int i, j;
    int s, t;    //子矩陣左上方元素行列下標
    int row, col;  //row用於存放階梯形矩陣非零行行數
    float swap;
    mc1 *head, *tail, *psnew;

    for (i=0; i<M; i++)    //輸入並初始化增廣矩陣
    {
        printf("請輸入增廣矩陣第%d行\n", i+1);
        scanf("%s", str);
        for (j=0; j<N; j++)
            a[i][j]=str[j]-48;
    }

    head=(mc1 *) malloc(sizeof(mc1));
    head->next=NULL;
    tail=head;
    s=t=1;  //子矩陣即為增廣矩陣本身，用增廣矩陣左上方元素行列標初始化s,t

    while((col=test(s, t, a))!=0)  //子矩陣不為零矩陣
    {
        if (s==M)  //增廣矩陣已化為階梯形矩陣
        {
            row=s;  //記錄非零行行數
            add(head, col, &tail);  //最後一個非零行主元列標放入maincol類型鏈表
            break;    //結束迴圈
        }
        else
        {
            j=s-1;
            for (i=s; i<M; i++)
            {
                if (fabs(a[j][col-1])<fabs(a[i][col-1]))    //列選主元
                    j=i;
            }

            if (s-1!=j)
            {
                for (i=col-1; i<N; i++)
                {
                    swap=a[j][i];
                    a[j][i]=a[s-1][i];    //列選主元
                    a[s-1][i]=swap;
                }
            }

            if (col==N)    //增廣矩陣已經化為階梯形矩陣
            {
                row=s;        //記錄非零行行數
                add(head, col, &tail);      //最後一個非零行主元列標放入maincol類型鏈表
                break;    //結束迴圈
            }

            for (i=s; i<M; i++)
                a[i][col-1]=-(a[i][col-1]/a[s-1][col-1]);

            for (i=col; i<N; i++)                        //消元
            {
                for (j=s; j<M; j++)
                    a[j][i]=a[j][col-1]*a[s-1][i]+a[j][i];
            }

            add(head, col, &tail);      //將消元後得到的主元列標col放入maincol類型鏈表
            s++;
            t=col+1;          //更新s,t,使s,t成為消元後得到的新的子矩陣的左上方元素行列標,為test函數操作新的子矩陣作準備
            continue;      //開始新一輪迴圈
        }
    }
    if (col==0)    //從迴圈控制條件退出迴圈
        row=s-1;  //此時增廣矩陣已成為階梯形矩陣,非零行函數就是s-1

    if (row==0)          //利用線性方程組解的判別準則判斷是否有解，有多少個解
    {
        printf("線性方程組有無窮多組解\n");  //增廣矩陣為零矩陣，無窮多組解
        printf("一般解:\n");
        for (i=1; i<N; i++)
            printf("x%d=t%d\n", i, i);  //輸出解
    }
    else
    {
        psnew=head->next;
        if (psnew->col==N)      //階梯形矩陣最後一主元在最後一列，無解
            printf("線性方程組無解\n");               
        else
        {
            convert(a, row, head);  //用convert函數將階梯形矩陣進一步化為簡化行階梯形矩陣
            if (row==N-1)    //非零行行數等於未知量個數，有唯一解
            {
                printf("線性方程組有唯一解:\n");
                for (i=1; i<=row; i++)                  //輸出唯一解
                    printf("x%d=%f\n", i, a[i-1][N-1]);
            }
            else    //非零行行數小於未知量個數，有無窮多組解
            {
                int *m;
                m=(int *) malloc((N-1-row)*sizeof(int));

                i=N-1-row;
                for (j=N-1; j>=1; j--)
                {
                    if (j!=psnew->col)
                    {
                        m[--i]=j;        //從所有未知量標號中篩選出自由未知量標號         
                        if (i==0)
                            break;
                    }
                    else
                    {
                        if (psnew->next!=NULL)
                            psnew=psnew->next;
                    }
                }
                printf("線性方程組有無窮多組解\n");
                printf("一般解:\n");
                i=row;
                for (psnew=head->next; psnew!=NULL; psnew=psnew->next)
                {
                    printf("x%d=%f", psnew->col, a[i-1][N-1]);    //輸出一般解
                    for (j=0; j<N-1-row; j++)
                    {
                        if(m[j]<psnew->col)
                        {
                            printf("-%dx%d", 0, m[j]);
                        }
                        else
                        {
                            printf("-%fx%d", a[i-1][m[j]-1], m[j]);
                        }
                    }
                    printf("\n");
                    i--;
                }
            }
        }
    }
}

int test(int s, int t, float a[][N])        //判斷增廣矩陣的s行至M行與t列至N列相交形成的子矩陣是否為零矩陣，若是返回0，若不是返回第一個不為零的列的列標
{
    int i, j;

    for (j=t-1; j<N; j++)
    {
        for (i=s-1; i<M; i++)
        {
            if (a[i][j]!=0)
                return (j+1);
        }
    }
    return (0);
}

void add(mc1 *head, int col, mc1 **tail)    //函數，用於建立一節點，其中儲存有主元col列標，然後按遞減順序將其插入maincol類型的鏈表
{
    mc1* psnew;

    psnew=(mc1 *) malloc(sizeof(mc1));
    psnew->col=col;

    if(head->next==NULL)
    {
        psnew->next=NULL;
        head->next=psnew;
        *tail=psnew;
    }
    else
    {
        psnew->next=head->next;
        head->next=psnew;
    }
}

void convert(float a[][N], int row, mc1 *head)    //函數，用於將階梯型矩陣轉化為簡化行階梯形矩陣
{
    mc1 *psnew, *pq;
    int i, j, k, m;

    psnew=head->next;
    for (i=row-1; i>=0; i--)
    {
        if (a[i][psnew->col-1]!=1)    //各非零行主元化為1
        {
          for (j=psnew->col; j<N; j++)
              a[i][j]=a[i][j]/a[i][psnew->col-1];
        }

        psnew=psnew->next;
    }

    psnew=head->next;          //向上消元把除第一個主元外其餘主元所在列中在主元上方的部分變為零
    for (i=row-1; i>=1; i--)
    {
        m=N-psnew->col-(row-i);  //擷取未知量標號1,2,--,N-1中位於i+1非零行主元列標號右邊的自由未知量標號個數
        for (j=i-1; j>=0; j--)
        {
            pq=head->next;    //pq指向存放最後一個非零行主元列標號的節點
            for (k=N; k>psnew->col; k--)
            {
                if (k!=pq->col)
                {
                  a[j][k-1]=-(a[i][k-1]*a[j][psnew->col-1])+a[j][k-1];  //從右向左作初等行變換直至i+1行主元所在列右邊的列位置，期間跳過i+2----row行主元所在的列
                  m--;
                  if (m==0)
                      break;
                }
                else
                {
                    if (pq->next!=psnew)
                        pq=pq->next;
                }
            }
        }
        psnew=psnew->next;  //遞進至上一行主元，為新一輪向上消元作準備
    }
}