考虑线性方程组

\[\mathrm{A}x=\mathrm{b}
\]

其中，\(\mathrm{A}=(a_{ij})_{n\times n}\)，\(\mathrm{b}=[b_1,b_2,\cdots,b_n]^{\mathrm{T}}\)。在线性代数的课程中，我们已经学习过Gauss消元法，具体操作是将矩阵A转化为“阶梯型”矩阵。为方便起见，本文仅仅讨论系数矩阵非奇异的方程组，此时，目标是将矩阵A转化为上三角矩阵，再执行回代过程，即可给出方程组的解。本文将给出在计算机上的具体操作及实例代码。

一、基本Gauss消去法

我们仅仅讨论对矩阵第一列的操作，剩余的操作可以以此类推，因而不再赘述。

在执行Gauss消去法时，我们将第一列对角元以下的元素全部变为零。记第一列消元操作后的增广矩阵为\([\mathrm{A}^{(1)},\mathrm{b}^{(1)}]\)，容易知道

\[[\mathrm{A}^{(1)},\mathrm{b}^{(1)}]=
\begin{bmatrix}
a_{11} & a_{22} & \cdots &a_{1n} & b_1 \\
0 & a_{22}^{(1)} &\cdots &a_{2n}^{(1)} & b_2^{(1)}\\
\vdots &\vdots & & \vdots &\vdots\\
0 & a_{n2}^{(1)} & &a_{nn}^{(1)} & b_n^{(1)}
\end{bmatrix}\]

其中

\[a_{ij}^{(1)}=a_{ij}-\frac{a_{i1}}{a_{11}}a_{1j}，j=2,\cdots ,n
\]
\[a_{i1}^{(1)}=0
\]
\[b_i^{(1)}=b_i-\frac{a_{i1}}{a_{11}}b_1
\]

观察到重复出现的结构\(\frac{a_{_{i1}}}{a_{_{11}}}\)，我们记它为\(l_{i1}\)，称为消元因子，并将它存储在原来\(a_{i1}\)的位置。在计算的过程中，先计算消元因子并存储在相应位置，再执行后续的算法。

对于后续部分的运算，在第k步，只要对矩阵\(A^{(k-1)}(k:n,k:n)\)执行相同操作即可。

二、列主元Gauss消去法

在执行Gauss消元法的过程中，如果\(a_{kk}^{(k-1)}\)相对于其他元素绝对值较小，则会产生较大的舍入误差，影响计算精度，为此，我们引入了列主元Gauss消去法，基于交换矩阵的行不影响线性方程组的解。

记执行完k-1步消元后的增广矩阵为\([\mathrm{A}^{(k-1)},\mathrm{b}^{(k-1)}]\)。考虑第k列对角元及其以下的部分。选择绝对值最大的元所在行，与当前行执行行交换，再进行Gauss消元法。

三、计算实例

用列主元Gauss消去法解以下线性方程组：

\[\left\{
\begin{array}{}
0.5x_1+1.1x_2+3.1x_3=6,\\
2x_1+4.5x_2+3.6x_3=0.02,\\
5x_1+0.96x_2+6.5x_3=0.96.
\end{array}
\right.\]

#include <iostream>

#include <math.h>

using namespace std;

int main()

{

    double A_Extended[3][4]={0.5,1.1,3.1,6,2,4.5,3.6,0.02,5,0.96,6.5,0.96};

    double X_solution[3];

    for (int i=0;i<=2;i++)

    {

        int n=i;

        for (int p=i+1;p<=2;p++)

        {

            if (fabs(A_Extended[p][i])>fabs(A_Extended[n][i]))

            {

                n=p;

            }

        }

        for (int p=i;p<=2+1;p++)

        {

            double k=A_Extended[n][p];

            A_Extended[n][p]=A_Extended[i][p];

            A_Extended[i][p]=k;

        }

        for (int p=i+1;p<=2;p++)

        {

            A_Extended[p][i]=-A_Extended[p][i]/A_Extended[i][i];

            for (int pco=i+1;pco<=2+1;pco++)

            {

                A_Extended[p][pco]=A_Extended[p][pco]+A_Extended[p][i]*A_Extended[i][pco];

            }

        }

    }

    X_solution[2]=A_Extended[2][3]/A_Extended[2][2];

    for (int i=1;i>=0;i--)

    {

        double sum=0;

        for (int k=2;k>i;k--)

        {

            sum=sum+A_Extended[i][k]*X_solution[k];

        }

        X_solution[i]=(A_Extended[i][3]-sum)/A_Extended[i][i];

    }

    cout<<X_solution[0]<<" "<<X_solution[1]<<" "<<X_solution[2]<<endl;

    return 0;

}

线性方程组的直接解法——Gauss消去法的相关教程结束。