[复习笔记]高斯消元法
概念
引用自百度百科：
数学上，高斯消元法（或译：高斯消去法），是线性代数规划中的一个算法，可用来为线性方程组求解。但其算法十分复杂，不常用于加减消元法，求出矩阵的秩，以及求出可逆方阵的逆矩阵。不过，如果有过百万条等式时，这个算法会十分省时。一些极大的方程组通常会用迭代法以及花式消元来解决。当用于一个矩阵时，高斯消元法会产生出一个“行梯阵式”。高斯消元法可以用在电脑中来解决数千条等式及未知数。亦有一些方法特地用来解决一些有特别排列的系数的方程组。
应用
由以上概念可知，高斯消元法可以用于求解线性方程组（即n元一次方程组），这也是高斯消元法在OI和ACM竞赛中的主要应用。当然也有求行列式等其他用法，但我太菜不会，学会后会更新
原理
首先回忆一下我们求解二元一次方程组的方法：代入消元与加减消元，其核心思想就是消元。实际上，求解n元一次方程组也是如此，都是通过加减与代入的方法实现消元
由于线性方程组的加减与倍乘都与矩阵的初等行变换类似，于是我们可以用矩阵来代替这一线性方程组。
我们可以将一个线性方程组的未知数省略，只留下系数和常数，并将其等效为一个增广矩阵（注：不会细说线代知识，无线代基础的选手可以放心食用）
以洛谷模板题)样例为例：
\(\left\{ \begin{matrix} x_1+3x_2+4x_3=5\\ x_1+4x_2+7x_3=4\\ 9x_1+3x_2+2x_3=2 \end{matrix} \right.\)
将其化为增广矩阵
\(\begin{bmatrix} 1&3&4&5\\ 1&4&7&3\\ 9&3&2&2 \end{bmatrix}\)
然后我们就可以通过矩阵的初等行变换（即方程的加减与倍乘）来实现消元，最终的目标是将其化为简化阶梯形矩阵（简单理解为除\(a_{ii}\)为1外其他系数都为0的矩阵）
即形如
\(\begin{bmatrix} 1&0&\cdots&0&b_{1}\\ 0&1&\cdots&0&b_{2}\\ \vdots&\vdots&\ddots&\vdots&\vdots\\ 0&0&\cdots&1&b_{n} \end{bmatrix}\)
易知此时的矩阵可以等价为
\(\left\{ \begin{matrix} x_1=b_1\\ x_2=b_2\\ \vdots\\ x_n=b_n \end{matrix} \right.\)
即原方程组的解集
算法过程
首先将矩阵化为阶梯形矩阵
也就是每一非零行都在每一零行之上，某一行的先导元素所在的列位于前一行先导元素的右边，某一先导元素所在列下方元素都是零的矩阵
具体为形如
\(\begin{bmatrix} a_{11}&a_{12}&\cdots&a_{1n}\\ 0&a_{22}&\cdots&a_{2n}\\ \vdots&\vdots&\ddots&\vdots\\ 0&0&\cdots&a_{mn} \end{bmatrix}\)
的矩阵
要实现这一过程的实现方法有很多，这里采用列主元消去法
对于第i个未知数，我们在还没有进行过消元的行（即第i行到第n行）中选取\(\left|a_i\right|\)最大那一个作为主元，将这一行交换到第i行并进行消元
你或许会发现，有时我们在第 i行到第n行找不到一个方程存在\(x_i\)即对于任意\(k \in [i,n]\)，都有\(a_{ki}=0\)。此时该线性方程组就无解或有无数组解。由于我们已经对第1行到第i-1行进行过消元，所以第i行到第n行都不存在\(x_1、x_2\cdots x_{i-1}\)。假如我们可以通过第i行到第n行这n-i+1行解出除\(x_i\)外剩下n-i个未知数，回代后第1行到第i-1行这i-1行中仍有i个未知数，这显然是无法得到唯一解的。
具体的消元过程为
- 选择\(\left|a_i\right|\)最大的一个为主元，将其所在行一行与当前行交换
- 将主元系数化为1
- 将第i+1到第n行减去第i行的\(a_{ki}\)倍，实现消去\(x_i\)
该部分代码如下
for(re int i=1;i<=n;++i)
{
re int p=i;
for(re int j=i+1;j<=n;++j)
if(fabs(mat[j][i])>fabs(mat[p][i]))
p=j;//选择绝对值最大的为主元
for(re int k=1;k<=n+1;++k) swap(mat[i][k],mat[p][k]);//交换两行
if(fabs(mat[i][i])<=eps){puts("No Solution");return 0;}
//找不到主元，无解或无数组解
re double tmp=mat[i][i];
for(re int j=i;j<=n+1;++j) mat[i][j]/=tmp;//化系数为1
for(re int j=i+1;j<=n;++j)
if(mat[j][i]!=0)
{
tmp=mat[j][i];
for(re int k=i;k<=n+1;++k)
mat[j][k]-=tmp*mat[i][k];
}
//实现消元
}
将以上过程完成后，我们就得到了阶梯形矩阵
接下来我们就可以进行回代。所谓回代，就是将已经求出来的未知数代回到上面的方程中。
由消元过程可知，最后一行\(a_1、a_2\cdots a_{n-1}=0\)。易知，该方程组有解，当且仅当最后一行，\(a_n=1\)，此时的\(b_n\)即为\(x_n\)的解，将其代入第1行到第n-1行就可以将其他方程的\(x_n\)消去，回代\(x_n\)后，我们便得到了\(x_{n-1}\)的解，重复以上过程就可以解出所有未知数
该部分代码如下
for(re int i=n;i>=1;--i)
for(re int j=i-1;j>=1;--j)
{
mat[j][n+1]-=mat[i][n+1]*mat[j][i];
mat[j][i]=0;
}
以上是高斯消元法的全过程
完整代码
#include <cmath>
#include <cstdio>
#include <algorithm>
using namespace std;
#define il inline
#define re register
const int N=110;
const double eps=1e-6;
namespace FastIO
{
char buf[1<<21],buf2[1<<21],*p1=buf,*p2=buf;
int p,p3=-1;
il int getc(){return p1==p2&&(p2=(p1=buf)+fread(buf,1,1<<21,stdin),p1==p2)?EOF:*p1++;}
il void flush(){fwrite(buf2,1,p3+1,stdout),p3=-1;}
#define isdigit(ch) (ch>=48&&ch<=57)
template <typename T>
il void read(T &x)
{
re bool f=0;x=0;re char ch=getc();
while(!isdigit(ch)) f|=(ch=='-'),ch=getc();
while(isdigit(ch)) x=x*10+(ch^48),ch=getc();
x=f?-x:x;
}
template <typename T>
il void print(T x)
{
if(p3>(1<<20)) flush();
if(x<0) buf2[++p3]=45,x=-x;
re int a[50]={};
do{a[++p]=x%10+48;}while(x/=10);
do{buf2[++p3]=a[p];} while (--p);
buf2[++p3]='\n';
}
}
using namespace FastIO;
int n;
double mat[N][N];
int main()
{
read(n);
for(re int i=1;i<=n;++i)
for(re int j=1;j<=n+1;++j)
read(mat[i][j]);
for(re int i=1;i<=n;++i)
{
re int p=i;
for(re int j=i+1;j<=n;++j)
if(fabs(mat[j][i])>fabs(mat[p][i]))
p=j;
for(re int k=1;k<=n+1;++k) swap(mat[i][k],mat[p][k]);
if(fabs(mat[i][i])<=eps){puts("No Solution");return 0;}
re double tmp=mat[i][i];
for(re int j=i;j<=n+1;++j) mat[i][j]/=tmp;
for(re int j=i+1;j<=n;++j)
if(mat[j][i]!=0)
{
tmp=mat[j][i];
for(re int k=i;k<=n+1;++k)
mat[j][k]-=tmp*mat[i][k];
}
}
for(re int i=n;i>=1;--i)
for(re int j=i-1;j>=1;--j)
{
mat[j][n+1]-=mat[i][n+1]*mat[j][i];
mat[j][i]=0;
}
for(re int i=1;i<=n;++i)
printf("%.2lf\n",mat[i][n+1]);
flush();return 0;
}
