线性方程组之高斯消元

一、形式

(egin{cases} A_{11}X_{1} + A_{12}X_{2} + ......+A_{1n}X_{n}=B_{1} \ A_{21}X_{1} + A_{22}X_{2} + ......+A_{2n}X_{n}=B_{2}\ ...... A_{n1}X_{1} + A_{n2}X_{2} + ......+A_{nn}X_{n}=B_{n} end{cases})

上列方程组：

(A=)(egin{bmatrix} A_{11} & cdots & A_{1n} \ vdots & vdots & vdots \ A_{n1} & cdots & A_{nn} end{bmatrix})

(B=)(egin{bmatrix} B_{1} \ vdots \ B_{n} end{bmatrix})

(X=)(egin{bmatrix} X_{1} \ vdots \ X_{n} end{bmatrix})

可以看作为：(AX=B)

二、初等变换（行）

1.将某行同乘或同除一个数（非0）
2.将某行加到另一行
3.将任意两行互换

值都不变（相信大家都知道）

三、步骤

大家可以发现，解线性方程其实相当于解小学方程，于是我们的步骤就是相当于模拟小学方程解法

1.消元

我们考虑对于第(i)个方程，我们消元消掉(X_{i})

于是我们需要找到每一个方程中，找到对于(X_i)的最大系数(maxn)以及其对应的方程，将它与第(i)个方程调换，这样能够保证最大系数的(X_i)处理次数最少，时间最少。

找到最大系数后，我们需要判断一下(maxn)是否等于0，如果等于0就代表没有解。但是因为使用(double)类型，要考虑下精度问题，所以判断条件为小于(eps)，

for(int i=1;i<=n;i++)
{
	int now=i;
	for(int j=i+1;j<=n;j++)
	{
		if(fabs(a[now][i])<fabs(a[j][i]))now=j;//寻找x[i]的最大系数
	}
	if(fabs(a[now][i])<eps)//判断，如果x[i]的系数小于eps，说明x[i]=0，判断为无解
	{
		puts("No Solution");
		return 0;
	}
	if(i!=now)swap(a[i],a[now]);//将最大系数x[i]所在的方程与第i方程交换 

我们由小学奥数可得，我们若想要把(X_i)消掉，最方便的方法是把所有方程的(X_i)的系数都化为(1)，再加减消元即可，在处理过程中需要注意方程中的其他变量也要随之变化

	dl div=a[i][i];//这就是maxn
	for(int j=i;j<=n+1;j++)a[i][j]/=div;//先处理以maxn为系数的x[i]所在的方程
	for(int j=i+1;j<=n;j++)//将其他方程都进行处理
	{
		div=a[j][i];//其他方程的x[i]的系数
		for(int k=i;k<=n+1;k++)
		{
			a[j][k]-=div*a[i][k];
		}		
	}
}

这样我们就可以一直循环，处理出最后的(X_n)

ans[n]=a[n][n+1]/a[n][n];

2.回代

回代也很好理解，每次都从(n)往回模拟方程代入，处理出每个(X_i)

for(int i=n-1;i>=1;i--)
{
	ans[i]=a[i][n+1];
	for(int j=i+1;j<=n;j++)ans[i]-=ans[j]*a[i][j];
}

整体代码

#include<bits/stdc++.h>
using namespace std;
typedef double dl;
const int N=110;
const dl eps=1e-9;
dl a[N][N],ans[N],x[N][N];
int n;

int main()
{
	scanf("%d",&n);
	for(int i=1;i<=n+1;i++)
	{
		for(int j=1;j<=n;j++)
		{
			scanf("%lf",&x[i][j]);
		}
	}
	for(int i=1;i<=n;i++)
	{
		for(int j=1;j<=n;j++)
		{
			a[i][j]=2*(x[1][j]-x[i+1][j]);
			a[i][n+1]+=x[1][j]*x[1][j]-x[i+1][j]*x[i+1][j];
		}
	}
	for(int i=1;i<=n;i++)
	{
		int now=i;
		for(int j=i+1;j<=n;j++)
		{
			if(fabs(a[now][i])<fabs(a[j][i]))now=j;
		}
		if(fabs(a[now][i])<eps)
		{
			puts("No Solution");
			return 0;
		}
		if(i!=now)swap(a[i],a[now]);
		dl div=a[i][i];
		for(int j=i;j<=n+1;j++)a[i][j]/=div;
		for(int j=i+1;j<=n;j++)
		{
			div=a[j][i];
			for(int k=i;k<=n+1;k++)
			{
				a[j][k]-=div*a[i][k];
			}		
		}
	}
	ans[n]=a[n][n+1]/a[n][n];
	for(int i=n-1;i>=1;i--)
	{
		ans[i]=a[i][n+1];
		for(int j=i+1;j<=n;j++)ans[i]-=ans[j]*a[i][j];
	}
	for(int i=1;i<=n;i++)printf("%.3lf ",ans[i]);
	return 0;
}

四、判断解

(A、有唯一解：)最后处理出来的(X_i)的系数不为0并且过程中得出系数不为0

(B、无解：)在处理过程中发现对于(X_i)的(maxn=0)

(C、无穷解：)最后处理出来的(X_i)的系数为0并且等式右边也为0

题目：
P3389 【模板】高斯消元法
P2455 [SDOI2006]线性方程组

五、总结

其实可以发现，高斯消元还是很好理解的，主要难的是前面的怎样推导线性方程组的过程，推导出来后就直接套模板就好了