【学习笔记】高斯消元法

P3389 【模板】高斯消元法

数论小白都能看懂的线性方程组及其解法

约旦消元法：

• 首先找到当前系数（i）的最大（主元）的柿子（为了尽可能减少爆精度），并调整到第i行
• 若当前最大项系数为0，No Solution
• 接下来进行对每个系数的非主元行的消元。形成俩全0三角形

#include<bits/stdc++.h>
using namespace std;
#define rep(i,a,b) for(int i=a;i<=b;++i)
#define dwn(i,a,b) for(int i=a;i>=b;--i)
template <typename T> inline void rd(T &x){x=0;char c=getchar();int f=0;while(!isdigit(c)){f|=c=='-';c=getchar();}while(isdigit(c)){x=x*10+(c^48);c=getchar();}x=f?-x:x;}

int n;
double a[110][110];

int main(){
	#ifdef WIN32 
	freopen("xiaoyuan.txt","r",stdin);
	#endif
	rd(n);
	rep(i,1,n)
		rep(j,1,n+1)
			scanf("%lf",&a[i][j]);
	rep(i,1,n){//枚举每一列的系数 
		int max=i;//假定这一列最大的系数是第i行的；
		rep(j,i+1,n)
			if(fabs(a[j][i])>fabs(a[max][i]))
				max=j;
		rep(k,1,n+1)
			swap(a[max][k],a[i][k]);//这样每次第i个未知数最大的系数都在第i行，可以保证俩个0三角形
		if(!a[i][i]){//最大值等于0则说明该列都为0，肯定无解 
			printf("No Solution");
			return 0;
		}
		rep(j,1,n){//每一行都减去一个数
			if(i!=j){//不是主元那一行
				double tmp=a[j][i]/a[i][i];//注意这里是浮点数！！！！ 
				rep(k,i,n+1)
				a[j][k]-=a[i][k]*tmp;
			}
		}
	}
	rep(i,1,n)
		printf("%.2lf
",a[i][n+1]/a[i][i]);
	return 0;
}
/*
3
1 3 4 5
1 4 7 3
9 3 2 2
*/
/*
-0.97
5.18
-2.39
*/