[笔记] 浅谈高斯消元

高斯消元。。。当初以为自己学会了。。。后来。。。

话说这个东西好像最早出现于《九章算术》诶（古代人就是强）

废话不说，进入正题。。。

前置知识

高斯消元法是解线性方程组的方法之一

首先，线性方程组了解一下：

　　可认为线性方程组就是一次方程组。如图：

如果存在常数c1,c2,c3,...,cn代替x1,x2,x3,...,xn,使上图的每个方程成立，则称（c1,c2,c3,...,cn）为方程组的一个解. 解方程就是求其全部解

显然对线性方程组进行以下三种变换,所得方程组解集不变(这三种变换被称为同解变换)

　　(1) 对换两个方程（换法变换）（其实就是交换两个方程的位置，上下交换) ;

　　(2) 用非零数 c 乘以某一个方程（倍法变换）;

　　(3) 把某一个方程的 k 倍加到另一个方程上去( 消法变换) .

的确显然吧...

再来了解一下矩阵(=_=)：

由 s*n 个数 a_ij ( i = 1, 2, … , s; j = 1, 2, … , n) 排成的矩形表(可理解为长是n，宽是s的二维数组。。)，称为 s 行 n 列矩阵 .a_ij称为矩阵的( i , j) 元 .通常用大写字母或 ( a_ij )_sn表示矩阵 .如果 s = n（长等于宽）称A是 n 阶方阵或 n 阶矩阵 .

全零矩阵记做O_sn,或O；

n阶单位矩阵记做I_nn,或I（单位矩阵是个方阵，从左上角到右下角的对角线（称为主对角线）上的元素均为1。除此以外全都为0。）

对于矩阵(a_ij)_sn、(b_ij)_ns，若a_ij=b_ji(1<=i<=s,1<=j<=n)则称B为A的转置矩阵，记做B=A^T

另：系数矩阵和增广矩阵

把线性方程组的系数a_ij列入矩阵，称为系数矩阵。

如图：

而在系数矩阵右侧加入一列常数项，称其为增广矩阵。

如图：

（好吧增广矩阵就是类似把未知数扣掉了再扔到矩阵里。。）

对矩阵的行(列)做线性方程组的三种同解变换（上面的），称为矩阵的初等变换。

高斯消元法

对方程组的增广矩阵做初等变换，将系数部分化为对角线，上三角形，或阶梯形

思路：找到一个不为零的系数，设其所在行为pos，列为i，然后拿这个数去消掉其它行，但同列的数；但是消元的时候要有一定顺序

步骤：

先把一个个方程扔到矩阵中，构建增广矩阵

然后从第i行（i从1到n，按顺序）开始，在i+1行到n中的第i列（上面的方程已经符合上三角性质了）选取一个不为零的系数，记其所在行为pos，

交换第pos行与第i行，然后pos行变成了第i行（注：后面的第i行为交换完后的，即交换前的pos行）

至于为什么要找最大的，这样有利于提高精度。。（学长说的。。）

于是拿第i行的第i项去消i+1行到n行的第i项

这样可以消成上三角矩阵；

然后从最后一行开始，向上消去不需要的系数（最后一行只有一个系数，可以直接计算答案）

代码（Luogu 模板）：

#include<cstdio>
#include<iostream>
#define R register int
using namespace std;
const double eps=1E-11;
inline int g() {
    R ret=0,fix=1; register char ch; while(!isdigit(ch=getchar())) fix=ch=='-'?-1:fix;
    do ret=ret*10+(ch^48); while(isdigit(ch=getchar())); return ret*fix;
}
int n;
double a[110][110];
inline bool ck0(double x) {return x<eps&&x>-eps;}
inline void Gauss() {
    for(R i=1,pos;i<=n;++i) {
        pos=i; for(R j=i+1;j<=n;++j) if(ck0(a[pos][i])||a[pos][i]<a[j][i]) pos=j;
        if(pos!=i) for(R j=1;j<=n+1;++j) swap(a[i][j],a[pos][j]);
        if(ck0(a[i][i])) {printf("No Solution
"); return ;}
        for(R j=i+1;j<=n;++j) {
            register double tmp=a[j][i]/a[i][i];
            for(R k=1;k<=n+1;++k) a[j][k]-=tmp*a[i][k];
        }
    } for(R i=n;i>=1;--i) {
        for(R j=n;j>i;--j) a[i][n+1]-=a[i][j]*a[j][n+1];//把这行系数代表的答案消掉
        a[i][n+1]/=a[i][i];
    } for(R i=1;i<=n;++i) printf("%.2lf
",a[i][n+1]);
}
signed main() {
    n=g(); for(R i=1;i<=n;++i) for(R j=1;j<=n+1;++j) a[i][j]=g();
    Gauss();
}

(例题)

Luogu P3389 【模板】高斯消元法

Luogu P4035 [JSOI2008]球形空间产生器

POJ1830

etc.

2019.05.10