高斯消元
这里介绍的是高斯-约旦消元法。
相对于传统的高斯消元，约旦消元法的精度更好、代码更简单，没有回带的过程。
约旦消元法大致思路如下：
1.选择一个尚未被选过的未知数作为主元，选择一个包含这个主元的方程。
2.将这个方程主元的系数化为1。
3.通过加减消元，消掉其它方程中的这个未知数。
4.重复以上步骤，直到把每一行都变成只有一项有系数。
注意，高斯消元可能会有精度问题，因此我们要取主元系数最大的式子去更新其他式子
例题一P3389 【模板】高斯消元法
在这里，无解和无穷解我们都输出No Solution
代码
(不够严谨的，无法具体判断无解还是无数解)
#include<bits/stdc++.h>
using namespace std;
const int maxn=105;
int n;
double a[maxn][maxn],ANS[maxn];
void gauss(){
for(int i=1;i<=n;i++){
int maxx=i;
for(int j=i+1;j<=n;j++)
if(fabs(a[maxx][i])<fabs(a[j][i]))maxx=j;
if(i!=maxx)swap(a[i],a[maxx]);
if(!a[i][i])continue;
for(int j=1;j<=n;j++){
if(j==i)continue;
double tmp=a[j][i]/a[i][i];
for(int k=i;k<=n+1;k++)
a[j][k]-=a[i][k]*tmp;
}
}
return ;
}
int main(){
cin>>n;
for(int i=1;i<=n;i++)
for(int j=1;j<=n+1;j++)
cin>>a[i][j];
gauss();
for(int i=1;i<=n;i++){
if(a[i][i]==0){
cout<<"No Solution\n";
return 0;
}
else
ANS[i]=a[i][n+1]/a[i][i];
}
for(int i=1;i<=n;i++)printf("%.2f\n",ANS[i]);
return 0;
}
例题二P2455 [SDOI2006]线性方程组
这一次，我们要具体判断无解与无穷解
题意：在高斯消元的基础上，判断是正常的还是“无解”还是“有无数解”。
简述一下高斯消元的思想：让第i个式子对应第i元，并用第i个式子消去其它式子上的第i元，最终使方程组呈现出类似一条对角线的样子，即第i 个式子只在第 i元和等式右边有值，其余都是0。
一般的高斯消元每次会选择该项系数绝对值最大的来对应第 i元。这样，如果发现有一元系数为0，就能说明方程没有正常的解了。
那么我们来研究，没有正常的解，究竟是“无解”还是“有无数解”。
一个\(naive\)的想法是，在完成消元后枚举每一项，如果发现系数是0则讨论等式右边：如果为0，则有无数解；如果不为0，则无解。
其中无解优先级更高，即如果一个方程有一项无解且有一项有无数解，那么判定无解。
交上去喜提90。
存在漏洞：
2
0 2 3
0 0 0
答案显然是0，但是我们的程序会说这是-1。
我们发现，我们可以根据第一个式子来求出第二元，但是程序会用它来算第一元，并且以后算第二元的时候不会再考虑该式，
因为我们认为它已经对应第一元了！
也就是说，我们的答案受到消元顺序影响。
注意到，我们的症结在于把还有用的式子放到了用不上它的地方，并且以后也不来看它是否可用。
那么，我们就来看它是否可用。
原本的高斯消元中，我们认为只有 i之后的式子是可用的，因为我们不用管无解还是无数解，系数为0直接判掉。
但这里，我们有可能会在系数为0之后继续做下去。这就是受到消元顺序影响的原因。
那么，可用的就不仅仅是之后的式子，还有之前系数为0的式子。
于是解决了。
代码
#include<bits/stdc++.h>
using namespace std;
const int maxn=105;
const double eps=1e-7;
int n;
double a[maxn][maxn],ANS[maxn];
void gauss(){
for(int i=1;i<=n;i++){
int maxx=i;
for(int j=1;j<=n;j++){//从i-n改为了1-n
if(fabs(a[j][j])>=eps&&j<i)continue;//增加了这句话
if(fabs(a[maxx][i])<fabs(a[j][i]))maxx=j;
}
if(i!=maxx)swap(a[i],a[maxx]);
if(!a[i][i])continue;
for(int j=1;j<=n;j++){
if(j==i)continue;
double tmp=a[j][i]/a[i][i];//这句话不能放进里面去，因为a[j][i]会修改
for(int k=i;k<=n+1;k++)
a[j][k]-=a[i][k]*tmp;
}
}
return ;
}
int main(){
cin>>n;
for(int i=1;i<=n;i++)
for(int j=1;j<=n+1;j++)
cin>>a[i][j];
gauss();
bool flag=false;
for(int i=1;i<=n;i++)if(fabs(a[i][i])<=eps&&fabs(a[i][n+1])>eps)flag=true;
if(flag){cout<<-1<<endl;return 0;}//先判断无解
flag=false;
for(int i=1;i<=n;i++)if(fabs(a[i][i])<=eps&&fabs(a[i][n+1])<=eps)flag=true;
if(flag){cout<<0<<endl;return 0;}//再判断无数解,并注意误差
for(int i=1;i<=n;i++)
ANS[i]=a[i][n+1]/a[i][i];
for(int i=1;i<=n;i++)printf("x%d=%.2f\n",i,ANS[i]);
return 0;
}
例题三P2447 [SDOI2010] 外星千足虫
高斯消元解异或方程组裸题
乍一看，这好像并不好做，最暴力的方法是O(2^n)枚举
但是有了高斯消元，我们就可以不这么暴力了，而是可以借助高斯消元的思想来解了
举个例子：
现在我们有一个方程组长这样：
（当然，给出的例子是很简单的了）
但是我们可以借助这种思想：我们首先观察第一个方程：我们希望只有这一个方程中含有\(x_1\)（这也是高斯消元的思想），这样在其他变量已知的情况下我们就可以求出\(x_1\)了
那么根据异或的性质，我们用第一个方程逐个与含有\(x_1\)的方程异或，就可以消去所有的\(x_1\)！
于是做出来的结果就长这样：
接下来就是第二个方程：我们要把\(x_2\)消掉，于是我们的方法同上
结果长这样：
最后一个方程
这就是高斯消元解异或方程组的过程及原理
但是这里会出现一个问题：“自由元”问题
按道理，给出n个方程组应该是可以解出每个答案的，但是由于异或运算特殊的性质，有时n个方程是难以解出每个值的
为什么？
举个例子：
不难发现，这组方程本身解就不唯一（显然，(1,0,1)与(0,1,0)都是该方程组的解）
（如果去探索本质的原因：上面的异或方程在逻辑上等同于以下两个方程：\(x_1==x_3\)且\(x_1!=x_2\)，那么考虑取值的方式，不难发现解不是唯一的）
那么如果用正常的方式来解，最后会变成这种形式：
这也能看出来，因为丢失了一个方程啊！
这时我们就称\(x_3\)是一个自由元，也就是他的值既可以为\(0\)又可以为\(1\)，此时该异或方程组有很多组解！
对于自由元问题，目前仅有的办法就是\(dfs\)，枚举每个自由元的状态然后回代求解
但是对于这道题，我们没有必要这么做
由于对于解不唯一的情况，他只要求输出一句话，所以我们只需在解不唯一的情况下（也就是解了全部方程仍不能得到唯一解时）特判即可
接下来我们就可以研究如何输出最小方程数了
我们直接像上面一样一个一个元解，不过我们现在不需要找到最大的元，因为这题不会出现精度问题
因此我们找到第一个式子，我们没有使用过，并且该元的系数不为0，记录下该位置
如果我们找不到，我们就返回解不唯一
否则我们用该式子更新所有式子
另外，如果我们使用常规的写法的时间复杂度是\(O(n^2m)\)的，会超时
因此我们可以使用\(bitset\)进行优化，会给复杂度除以32，所以最终的时间复杂度是\(O(\frac{n^2m}{32})\)
代码
#include<bits/stdc++.h>
using namespace std;
const int maxn=1e3+5;
int n,m,maxx=0;
bitset<maxn>a[2*maxn];
string x;
char y;
bool Gauss(){
for(int i=1;i<=n;i++){
for(int j=1;j<=m;j++){
if(a[j][j]&&j<i||!a[j][i])continue;
if(i!=j)
swap(a[i],a[j]);
maxx=max(maxx,j);
break;
}
if(!a[i][i])return false;
for(int j=1;j<=m;j++){
if(i==j)continue;
if(a[j][i])
for(int k=i;k<=n+1;k++)
a[j][k]=(a[j][k]^a[i][k]);
}
}
return true;
}
int main(){
cin>>n>>m;
for(int i=1;i<=m;i++){
cin>>x>>y;
for(int j=1;j<=n;j++)a[i][j]=x[j-1]-'0';
a[i][n+1]=y-'0';
}
if(!Gauss()){
cout<<"Cannot Determine\n";
return 0;
}
cout<<maxx<<endl;
for(int i=1;i<=n;i++){
if(a[i][n+1])
printf("?y7M#\n");
else printf("Earth\n");
}
return 0;
}