高斯消元法&&列主消元法

写在前面

最近在学数值分析，正好学到求解线性方程组。就自己动手简单实现了一下。关于本算法的原理可以在《数值分析》第3版（马东升 董宁编），对应于该书的P145页，详细讲解了公式。因本人时间有限，暂时不详细编辑公式，等空闲了再来重新补充。简要的说明下该算法的应用吧，高斯消元法在线性代数那门课肯定学过了，对于一般简单的3、4阶线性方程组，还可以进行纸上的笔运算，但是当数目过多，算起来就比较吃力了，所以借助于计算机实现。如果很多的数学软件里都有矩阵的求解的函数，如mathematic和matlab，因本人对这些软件不熟悉，所以也没去去尝试使用。对于非奇异矩阵并且主对角线上的元素不为0。就可以用高斯消元法来求解。

高斯消元法

话不多说，直接贴代码

#include <bits/stdc++.h>
using namespace std;
const int n = 3;
void gaussin(double a[n][n], double b[n])
{
    for (int i = 0; i < n; i++)
    {
        if (a[i][i] == 0)
        {
            cout << "can't use gaussin meathod" << endl;
            return;
        }
    }
    int i, j, k;
    double c[n];
    for (k = 0; k < n - 1; k++)
    {
        for (i = k + 1; i < n; i++)
            c[i] = a[i][k] / a[k][k];
        for (i = k + 1; i < n; i++)
        {
            for (j = 0; j < n; j++)
            {
                a[i][j] = a[i][j] - c[i] * a[k][j];
            }
            b[i]=b[i]-c[i]*b[k];
        }
    }
    double x[n];
    for(i=n-1;i>=0;i--)
    {
        double sum = 0;
        for(j=i+1;j<n;j++)
        {
            sum+=a[i][j]*x[j];
        }
        x[i]=(b[i]-sum)/a[i][i];
    }
    cout<<"the solution of the equation is : "<<endl;
    cout<<endl;
    for(int i =0;i<n;i++)
    cout<<"x"<<i+1<<" = "<<x[i]<<endl;
}
int main()
{
    double a[3][3]={12,-3,3,-18,3,-1,1,1,1};
    double b[3]={15,-15,6};
    gaussin(a,b);
    return 0;
}    

列主消元法

高斯消元法,但是，在消元过程中，无法使主元素a(ii)≠0，但是很小时，用其做除数，会导致其他元素数量级的严重增长，舍入误差的扩展，最后导致计算结果不可靠。所以这次采用列主元素消去法来进行，思想就是将有小数的那行与该列中数最大的那行进行交换。 还是写帖代码吧

#include <bits/stdc++.h>
using namespace std;
const int n = 3;
void gaussin(double a[n][n], double b[n])
{
    int i, j, k;
    int col, row;
    double c[n];
    for (k = 0; k < n - 1; k++)
    {
        double ave = 0;
        for (i = k; i < n; i++)
        {
            if (fabs(a[i][k]) > ave)
            {
                ave = fabs(a[i][k]);
                row = i;
                col = k;
            }
        }
        if (a[row][row] == 0)
            cout << "can't solve" << endl;
        if (k != row)
        {
            for (i = 0; i < n; i++)
            {
                swap(a[row][i], a[k][i]);
            }
            swap(b[k], b[row]);
        }
        for (i = k + 1; i < n; i++)
            c[i] = a[i][k] / a[k][k];
        for (i = k + 1; i < n; i++)
        {
            for (j = 0; j < n; j++)
            {
                a[i][j] = a[i][j] - c[i] * a[k][j];
            }
            b[i] = b[i] - c[i] * b[k];
        }
    }
    double x[n];
    for (i = n - 1; i >= 0; i--)
    {
        double sum = 0;
        for (j = i + 1; j < n; j++)
        {
            sum += a[i][j] * x[j];
        }
        x[i] = (b[i] - sum) / a[i][i];
    }
    cout << "the solution of the equation is : " << endl;
    cout << endl;
    for (int i = 0; i < n; i++)
        cout << "x" << i + 1 << " = " << x[i] << endl;
}
int main()
{
    double a[3][3] = {12, -3, 3, -18, 3, -1, 1, 1, 1};
    double b[3] = {15, -15, 6};
    gaussin(a, b);
    return 0;
}