高斯消元法（列选主元法）
唯一解：判定存在性并求值
- a数组存增广矩阵（第 n+1列为常数列）
- mr为每次所选取的最大系数的行标号，然后交换第 i行和 mr行
- 若存在唯一解，则系数矩阵会被化简为 n×n的单位矩阵，相应的常数列就是 xi的解啦！
- 调用：Gauss();
- 返回值： 0表示无唯一解， 1表示有唯一解
#include "bits/stdc++.h"
#define hhh printf("hhh\n")
#define see(x) (cerr<<(#x)<<'='<<(x)<<endl)
using namespace std;
typedef long long ll;
typedef pair<int,int> pr;
inline int read() {int x=0;char c=getchar();while(c<'0'||c>'9')c=getchar();while(c>='0'&&c<='9')x=x*10+c-'0',c=getchar();return x;}
const int maxn = 1e5+10;
const int inf = 0x3f3f3f3f;
const int mod = 1e9+7;
const double eps = 1e-7;
int n;
double a[101][102];
bool Gauss() {
for(int i=1; i<=n; ++i) {
int mr=i;
for(int j=i+1; j<=n; ++j) if(fabs(a[j][i])>fabs(a[mr][i])) mr=j;
swap(a[i],a[mr]);
if(fabs(a[i][i])<eps) return 0;
for(int j=1; j<=n; ++j) {
if(j==i) continue;
double d=a[j][i]/a[i][i];
for(int k=i; k<=n+1; ++k) a[j][k]-=d*a[i][k];
}
for(int k=i+1; k<=n+1; ++k) a[i][k]/=a[i][i]; a[i][i]=1;
}
return 1;
}
int main() {
//ios::sync_with_stdio(false); cin.tie(0);
n=read();
for(int i=1; i<=n; ++i)
for(int j=1; j<=n+1; ++j) scanf("%lf", &a[i][j]);
if(Gauss()) for(int i=1; i<=n; ++i) printf("%.2f\n", a[i][n+1]);
else printf("No Solution\n");
}
无穷解：判定解的情况
- 要先判定无解，再判定无穷解（但是要在化简过程中先记录）
- 判定无解要在化简结束后，并且要检查所有行
#include "bits/stdc++.h"
#define hhh printf("hhh\n")
#define see(x) (cerr<<(#x)<<'='<<(x)<<endl)
using namespace std;
typedef long long ll;
typedef pair<int,int> pr;
inline int read() {int x=0;char c=getchar();while(c<'0'||c>'9')c=getchar();while(c>='0'&&c<='9')x=x*10+c-'0',c=getchar();return x;}
const int maxn = 1e5+10;
const int inf = 0x3f3f3f3f;
const int mod = 1e9+7;
const double eps = 1e-7;
int n;
double a[101][102];
int Gauss() {
bool Inf=0;
for(int i=1; i<=n; ++i) {
int mr=i;
for(int j=i+1; j<=n; ++j) if(fabs(a[j][i])>fabs(a[mr][i])) mr=j;
swap(a[i],a[mr]);
if(fabs(a[i][i])<eps) { Inf=1; continue; }
for(int j=1; j<=n; ++j) {
if(j==i) continue;
double d=a[j][i]/a[i][i];
for(int k=i; k<=n+1; ++k) a[j][k]-=d*a[i][k];
}
for(int k=i+1; k<=n+1; ++k) a[i][k]/=a[i][i]; a[i][i]=1;
}
for(int i=1; i<=n; ++i) {
int f=1;
for(int j=1; j<=n; ++j) if(fabs(a[i][j])>=eps) { f=0; break; }
if(f&&fabs(a[i][n+1])) return 0;
}
if(Inf) return 2;
return 1;
}
int main() {
//ios::sync_with_stdio(false); cin.tie(0);
n=read();
for(int i=1; i<=n; ++i)
for(int j=1; j<=n+1; ++j) scanf("%lf", &a[i][j]);
int ans=Gauss();
if(ans==1) for(int i=1; i<=n; ++i) printf("x%d=%.2f\n", i, a[i][n+1]);
else printf("%d\n", ans?0:-1);
}
高级版（可求出自由变量数目）（还没用过这个板子，但感觉可以用上面的写法稍微改改就好了）
a数组中常系数保存在a[i][var]中；
唯一解保存在x数组中；
自由变量保存在free_x数组中；
const double eps = 1e-9;
double a[maxn][maxn];//每个方程中每个未知数的系数
double x[maxn];//解集
bool free_x[maxn];
int sgn(double x) {
return (x>eps)-(x<-eps);
}
int gauss(int equ, int var) {
int i,j,k;
int max_r; // 当前这列绝对值最大的行.
int col; // 当前处理的列.
double temp;
int free_x_num;
int free_index;
// 转换为阶梯阵.
col=0; // 当前处理的列.
memset(free_x,1,sizeof(free_x));
for(k=0; k<equ&&col<var; k++,col++)
{
max_r=k;
for(i=k+1; i<equ; i++) {
if(sgn(fabs(a[i][col])-fabs(a[max_r][col]))>0) max_r=i;
}
if(max_r!=k){ // 与第k行交换.
for(j=k; j<=var; j++) swap(a[k][j],a[max_r][j]);
}
if(sgn(a[k][col])==0) { k--; continue; }// 说明该col列第k行以下全是0了，则处理当前行的下一列.
for(i=k+1; i<equ; i++) { // 枚举要删去的行.
if(sgn(a[i][col])!=0) {
temp=a[i][col]/a[k][col];
for(j=col; j<=var; j++) {
a[i][j]=a[i][j]-a[k][j]*temp;
}
}
}
}
for(i=k;i<equ;i++) if(sgn(a[i][col])!=0) return 0;
if(k<var) {
for(i=k-1; i>=0; i--) {
free_x_num=0;
for(j=0; j<var; j++) {
if(sgn(a[i][j])!=0&&free_x[j]) free_x_num++, free_index=j;
}
if(free_x_num>1) continue;
temp=a[i][var];
for(j=0; j<var; j++) {
if(sgn(a[i][j])!=0&&j!=free_index) temp-=a[i][j]*x[j];
}
x[free_index]=temp/a[i][free_index];
free_x[free_index]=0;
}
return var-k;
}
for(i=var-1; i>=0; i--) {
temp=a[i][var];
for(j=i+1; j<var; j++) {
if(sgn(a[i][j])!=0) temp-=a[i][j]*x[j];
}
x[i]=temp/a[i][i];
}
return 1;
}
京公网安备 11010502036488号