高斯消元法

高斯消元

这里只是概述，更详细的内容看一本通金牌导航吧，蒟蒻能力有限qwq

高斯消元，用于求解n元一次方程组
这里介绍高斯约旦消元法（反正我看的题解是这么叫的）
该方法的思路是：通过一系列消元操作，使得每一变量仅出现在一个方程中
简要流程如下：
(1)对于每一个n元一次方程，将形如∑a[i]x[i]=b的方程转化成行列式形式
(2)枚举每一行，假设我们枚举到第i行，我们希望x[i]最终仅在第i个方程中出现，称x[i]为该行的主元。如果第i行的主元系数不是第i~n个方程中最大的，则找到那个系数最大的方程，将行列式中对应的两行互换，使得第i行的主元系数最大化，以减小精度误差
(3)设主元系数为K，将第i行每一项系数除以K，目的是将主元系数变为1
(4)枚举除第I行以外的每一行，令该行x[i]的系数为K,将该行与第i行放大K倍得到的式子相减，消去该行的x[i]项
（5）重复上述步骤，直至n元一次不等式组化为n个一元一次方程，然后直接求解就行了

看了代码就会更明白一些
洛谷P4035 [JSOI2008]球形空间产生器

#include<bits/stdc++.h>
using namespace std;
#define il inline
#define re register
il int read(){
	int s=0,w=1;char c=getchar();
	while(c<'0'||c>'9'){ if(c=='-') w=-1;c=getchar();}
	while(c>='0'&&c<='9'){ s=(s<<1)+(s<<3)+c-'0';c=getchar();}
	return s*w;
}
const int N=20;
const double eps=1e-8;
int n,m;//m=n+1
double f[N][N];
il void Gauss(){
	for(re int i=1;i<=n;i++){
		int l=i;double K=0;
		for(re int j=i+1;j<=n;j++)
			if(fabs(f[j][i])-fabs(f[l][i])>eps) l=j;
		if(l!=i) for(re int j=i;j<=m;j++) swap(f[i][j],f[l][j]);
		K=f[i][i];
		for(re int j=i;j<=m;j++) f[i][j]/=K;
		for(re int j=1;j<=n;j++){
			if(j==i) continue;
			K=f[j][i];
			for(re int k=i;k<=m;k++) f[j][k]-=f[i][k]*K;
		}
	}
}
int main()
{
	n=read(),m=n+1;
	for(re int i=0;i<=n;i++){
		for(re int j=1;j<=n;j++) scanf("%lf",&f[i][j]);
	}
	for(re int i=1;i<=n;i++){
		double D=0;
		for(re int j=1;j<=n;j++){
			D+=f[i][j]*f[i][j]-f[0][j]*f[0][j];
			f[i][j]=2*(f[i][j]-f[0][j]);
		}
		f[i][m]=D;
	}
	Gauss();
	for(re int i=1;i<n;i++) printf("%.3lf ",f[i][m]/f[i][i]);
	printf("%.3lf",f[n][m]/f[n][n]);
	
}

end