高斯消元

算法分析

如果给定一个形如以下式子的多元方程式

\begin{cases} 2x+y-z=8\\ -3x-y+2z=11\\-2x+y+2z=-1\end{cases}⎩⎪⎨⎪⎧​2x+y−z=8−3x−y+2z=11−2x+y+2z=−1​

我们要首先提出各项的系数

• 因为我们知道，高斯消元其实只跟系数有关

我们可以写成以下的矩阵形式

其中左边是各项的系数，分隔线之后的是等式右边的常数列

\begin{bmatrix}2&1&-1&|&8 \\ -3&-1&2&|&-11\\-2&1&2&|&-3\end{bmatrix}⎣⎡​2−3−2​1−11​−122​∣∣∣​8−11−3​⎦⎤​

\quad\qquad\Downarrow⇓ 经过r行和第i行交换

\begin{bmatrix}-3&-1&2&|&-11 \\ 2&1&-1&|&8\\-2&1&2&|&-3\end{bmatrix}⎣⎡​−32−2​−111​2−12​∣∣∣​−118−3​⎦⎤​

此时就可以开始加减消元了

首先我们要用第ii个方程的来消去第kk个方程的第kk列，那么第kk行的所有元素 A[k][j]A[k][j]都应该减去A[i][j]A[i][j]的A[k][i]/A[i][i]A[k][i]/A[i][i]倍。

（我们这里约定，A[a][b]指的是第aa行第bb列的系数）

我们用实例讲解一下：

• 这是一个已经经过第一步处理的矩阵（我承认就是从上面copy过来的）：

\begin{bmatrix}-3&-1&2&|&-11 \\ 2&1&-1&|&8\\-2&1&2&|&-3\end{bmatrix}⎣⎡​−32−2​−111​2−12​∣∣∣​−118−3​⎦⎤​

\quad\qquad\Downarrow⇓ 经过加减消元

\begin{bmatrix}-3&-1&2&|&-11 \\&&&|\\ &\dfrac{1}{3}&\dfrac{1}{3}&|&\dfrac{2}{3}\\&&&|\\&\dfrac{5}{3}&\dfrac{2}{3}&|&\dfrac{13}{3}\end{bmatrix}⎣⎢⎢⎢⎢⎢⎢⎡​−3​−131​35​​231​32​​∣∣∣∣∣​−1132​313​​⎦⎥⎥⎥⎥⎥⎥⎤​

（这就是加减消元的第一步）

还不理解？没关系（我再解释一下）

• 请看这里

  1、第二行相当于是L_1*\dfrac{2}{3}+L_2L1​∗32​+L2​

  2、同理，第三行相当于是L_1*(-\dfrac{2}{3})+L_2L1​∗(−32​)+L2​

  （这里的L_kLk​指的是矩阵中的第kk行）

现在应该明白了吧

接下来要做什么呢？ \qquad——代入！

根据我们的算法，可以直接从最后一行推出X_nXn​的值，然后倒数第二行可以推出X_{n-1}Xn−1​的值，· · ·以此类推。

最后我们可以求出所有XX的唯一解

然后？. . . . . . 然后就没有然后了呀

至此，对于高斯消元的理论分析已全部结束，接下来给代码。

//主元： 
double gauss()
{
    double x=1;
    for(int a=1;a<=n;a++)
    {
        for(int b=a;b<=n;b++)
            if(fabs(z[b][a])>1e-8)
            {
                if(b==a)
                    break;
                x=-x;
                for(int c=1;c<=n;c++)
                    swap(z[b][c],z[a][c]);
                break;
            }
        if(fabs(z[a][a])<=1e-8)
            return 0.0;
        for(int b=a+1;b<=n;b++)
        {
            double z[b][a]/z[a][a];
            for(int c=1;c<=n;c++)
                z[b][c]=z[b][c]-z[a][c]*k;
        }
    }
    for(int a=1;a<=n;a++)
        x=x*z[a][a];
    return x;
}



//辗转相消：
int gauss()
{
    int x=1;
    for(int a=1;a<=n;a++)
    {
        for(int b=a+1;b<=n;b++)
        {
            while(z[b][a]!=0)//z是double类型 
            {
                int k=z[a][a]/z[b][a];
                for(int c=1;c<=n;c++)
                    z[a][c]=z[a][c]-z[b][c]*k;
                for(int c=1;c<=n;c++)
                    swap(z[a][c],z[b][c]);
                x=-x;
            }
        }
    }
    for(int a=1;a<=n;a++)
        x=x*z[a][a];
    return x;
}