【数论】高斯消元

高斯消元

定义

数学上，高斯消元法（或译：高斯消去法），是线性代数规划中的一个算法，可用来为线性方程组求解。但其算法十分复杂，不常用于加减消元法，求出矩阵的秩，以及求出可逆方阵的逆矩阵。不过，如果有过百万条等式时，这个算法会十分省时。一些极大的方程组通常会用迭代法以及花式消元来解决。当用于一个矩阵时，高斯消元法会产生出一个“行梯阵式”。高斯消元法可以用在电脑中来解决数千条等式及未知数。亦有一些方法特地用来解决一些有特别排列的系数的方程组。

——Baidu

高斯消元法（$Gauss-Jordan elimination$）是求解线性方程组的经典算法。并且高斯消元还可以计算行列式，求矩阵的逆。

消元法

消元法主要为加减消元，同时加减消元也是一个比较常用的方法，如：

$$\begin{equation} \left\{ \begin{aligned} 4x+y=10 ------(1)\\ x-y=5 ------ (2) \end{aligned} \right. \end{equation}$$

(1)+(2)式可得：

$$5x=15$$

解得：

$$x=3$$

将$x=3$带回(1)和(2)任一式子可得：

$$y=-2$$

这就是消元法，通过使两个方程中同一未知数的系数相加减来消去未知数，完成求解的目的。

高斯消元法

对于一个线性方程组我们可以将其化为一个矩阵，比如一个$m$个$n$元一次方程组，我们可以将其化为一个$m \times (n+1)$的一个矩阵，矩阵的每行的前$n$个元素表示为每个未知数的系数，第$n+1$个元素表示为该方程所等于的常数，这个矩阵就是增广矩阵。
如:
对于一个三元一次方程组：

$$x_1+x_2+x_3=6 \\ 2x_1-x_2+x_3=1 \\ x_1+2x_2-3x_3=1$$

其增广矩阵为：

$$\begin{pmatrix}1&1&1&|&6\\2&-1&1&|&1\\1&2&-3&|&1\end{pmatrix}$$

接下来会经过以下变化：

$$\begin{pmatrix}1&1&1&|&6\\2&-1&1&|&1\\1&2&-3&|&1\end{pmatrix}$$

$$\begin{pmatrix}1&1&1&|&6\\0&-3&-1&|&-11\\1&2&-3&|&1\end{pmatrix}$$

$$\begin{pmatrix}1&1&1&|&6\\0&-3&-1&|&-11\\0&1&-4&|&-5\end{pmatrix}$$

$$\begin{pmatrix}1&1&1&|&6\\0&-3&-1&|&-11\\0&0&-13&|&-26\end{pmatrix}$$

根据最后一个式子即可得：

$$x_3=2$$

接着回带：

$$\begin{pmatrix}1&0&0&|&1\\0&1&0&|&3\\0&0&1&|&2\end{pmatrix}$$

以上的几种变换叫做初等行变换，后化简后的叫做阶梯型矩阵，最后化简后的矩阵为简化阶梯型矩阵。
同时如果化简后出现这种情况：

$$x_1=6-x_2\\ x_3=1$$

这时$x_2$取任何一个值都有一个$x_1$与之对应，在上述的方程中，像$x_1,x_3$称为主元，$x_2$称为自由元。

引论

根据以上分析，高斯消元后会有三种情况：

• 若存在系数全为零，然而常数不为零时，方程无解；
• 若系数不全为零的行数和未知数的个数相等，则方程有唯一解；
• 若系数不全为零的行数小于未知数的个数，则方程有无数个解。

代码实现

首先我们整理一下高斯消元的步骤：

• 增广矩阵初等行变换为最简阶梯矩阵；
• 还原线性方程组；
• 求解第一个变量；
• 补充自由未知量；
• 列表示方程组通解。

CODE：

#include<bits/stdc++.h>
#define eps 1e-8
using namespace std;
int n;
double c[110][110];
int main()
{
    scanf("%d",&n);
    for(int i=1;i<=n;i++)
        for(int j=1;j<=n+1;j++)
            scanf("%lf",&c[i][j]);
    for(int i=1;i<=n;i++)
    {
        for(int j=i;j<=n;j++)
        {
            if(fabs(c[j][i])>eps)swap(c[i],c[j]); // 找到x[i]系数不为零的方程
        }
        for(int j=1;j<=n;j++)
        {
        	if(i==j)continue;
            double rate=c[j][i]/c[i][i];
            for(int k=i;k<=n;k++)
            {
                c[j][k]-=c[i][k]*rate;
                if(fabs(c[j][k])<eps)c[j][k]=0;// 浮点型不可避免的精度误差
            }//  消元
            c[j][n+1]-=c[i][n+1]*rate;
            if(fabs(c[j][n+1])<eps)c[j][n+1]=0;// 改变常数
        }
    }
    for(int i=1;i<=n;i++)
    {
    	int x=0;
    	for(int j=1;j<=n;j++)
    	{
    		if(c[i][j]==0)x++;
		}
		if(x==n)
		{
			printf("No Solution");
			return 0;
		}// 判断是否无解
	}
    for(int i=1;i<=n;i++)printf("%0.2lf\n",c[i][n+1]/c[i][i]);
    return 0;
}

附上模板题P3389 【模板】高斯消元法。

练习题：P4035 [JSOI2008]球形空间产生器。

完结撒花（>-<)