高斯消元法

一、基本描述

1.高斯消元

主要用来求解线性方程组，也可以求解矩阵的秩，矩阵的逆。在ACM中是一个有力的数学武器。

时间复杂度是O（n³），主要与方程组的个数，未知数的个数有关。

2.什么是线性方程组？

有多个未知数，并且每个未知数的次数均为一次，这样多个未知数组成的方程组为线性方程组。

二、算法实现

       高斯消元法就是手解方程组。

基本原理就是：消元（因为线性方程组次数为1次，不需要考虑降次），然后回代，在得到kx=b这类的式子后，就可以显而易见得到x的解了。

 

举例说明高斯消元法。

假设有如下的方程组：

写成矩阵形式就是： AX=B ，其中：

且 ，现对矩阵 A 作初等变换，同时矩阵 B 也作同样的初等变换，则当 A 化为单位矩阵的时候，有：

显而易见，我们得到了方程组的解：  。（Ps.这里的矩阵带T代表是列向量，即x1=1,x2=2,x3=4）

所以， 我们要以一定的策略， 对 A 和 B 施以一系列的初等变换， 当 A 化为单位矩阵的时候，B 就为方程组的解。

三、高斯消元法实现

网上很多模板还强调解出矩阵的上三角形，但是对于一般的题目来说，并没有必要。学会怎么解就好了。

直接例题上板子。

1.bzoj1013 [JSOI2008]球形空间产生器sphere

【题目大意】

给出n维空间中给出n+1个点的坐标，求出球心坐标。

【思路】

令球心坐标为x1,x2...xn,假设当前第i个点坐标为a1,a2...,an，第i+1个点坐标为b1,b2...,bn，则由半径相等可得：

(a1-x1)^2+(a2-x2)^2+...+(an-xn)^2=(b1-x1)^2+(b2-x2)^2+...+(bn-xn)^2

化简可得:

2(a1-b1)x1+2(a2-b2)x2+...+2(an-bn)xn=(a1^2+a2^2+...+an^2-b1^2-b2^2-...-b3^2)

如此可得到n个一元n次方程组，用最简单的高斯消元搞一搞就好了。

#include<iostream>
#include<cstdio>
#include<algorithm>
#include<cmath>
#define maxn 25
using namespace std;
const double eps=1e-6; 
double f[maxn],a[maxn][maxn];
int n;
void Gauss(){
    int i,j,now=1;//now表示处理到第几行
    for (i=1;i<=n;i++){  //化简为行阶梯型 
        for (j=now;j<=n;j++){
            if (fabs(a[j][i])>eps) break;  //当前未知量系数不为0
        }
        if (j>n) continue;
        if (j!=now){
            for (int k=1;k<=n+1;k++) swap(a[j][k],a[now][k]);
        }
        double t=a[now][i];
        for (int k=1;k<=n+1;k++) a[now][k]/=t;
        for (int k=1;k<=n;k++){
            if (k!=now){
                t=a[k][i];
                for (int l=1;l<=n+1;l++){
                    a[k][l]-=t*a[now][l];
                }
            }
        }
        now++;
    }
    return ;
}
int main(){
    double x;
    scanf("%d",&n);
    for (int i=1;i<=n;i++) scanf("%lf",&f[i]);
    for (int i=1;i<=n;i++){
        for (int j=1;j<=n;j++){
            scanf("%lf",&x);
            a[i][j]=2*(x-f[j]);
            a[i][n+1]+=x*x-f[j]*f[j]; 
        }
    }
    Gauss();
    for (int i=1;i<n;i++) printf("%.3lf ",a[i][n+1]);
    printf("%.3lf
",a[n][n+1]);
    return 0;
}