最近遇到了几道高斯消元的题，虽然学了线性代数知道求解方法，但无法通过编程实现，相信还有少部分人也遇到了这种问题，所以写下这篇博文
数学上，高斯消元法（或译：高斯消去法），是线性代数规划中的一个算法，可用来为线性方程组求解。但其算法十分复杂，不常用于加减消元法，求出矩阵的秩，以及求出可逆方阵的逆矩阵。不过，如果有过百万条等式时，这个算法会十分省时。一些极大的方程组通常会用迭代法以及花式消元来解决。当用于一个矩阵时，高斯消元法会产生出一个“行梯阵式”。高斯消元法可以用在电脑中来解决数千条等式及未知数。亦有一些方法特地用来解决一些有特别排列的系数的方程组。
简单的说，通过高斯消元我们可以方便的解决多元一次方程组。
2、思路
其实就是通过编程来模拟现实生活中解决线性代数问题。依次让某一列从第几行开始变换成0，得到上三角形矩阵，就可以判断出答案是无解、无穷多解还是有唯一解了。
大家可以看我待会前两题的代码理解一下，第一题代码有注释，最好自己动手模拟一下。
实在不懂的话，这里建议y总的acwing基础课套餐来一个或者去学一下线性代数。
不难，用心学就会了。
废话不多说，直接上模板和题目
3、题目
第一题： 题目来源：acwing里面的–高斯消元解线性方程组
#include<bits/stdc++.h>
using namespace std;
const int N=105;
const double eps=1e-6;
int n;
double a[N][N];
int guass(){
int c,r; //c为当前考虑的列数，即当前将第c列变成系数为1的行列式。 r为在考虑当前第c列的情况下，初始行数
for(c=0,r=0;c<n;c++){
int t=r; //t为在考虑当前第c列的情况下，系数最大的行数
for(int i=r;i<n;i++){
if(fabs(a[i][c])>fabs(a[t][c])){ //在当前考虑第c列的情况下，如果第i行第c列的数>第t行第c列的数 ，则t=i，用来记录最大的数是哪一行
t=i;
}
}
if(fabs(a[t][c])<eps) continue; //如果在考虑当前第c列的情况下，最大的那一行的第c列的系数为0，那么继续，直接结束该次循环，去之星下一列的重复操作
for(int i=c;i<=n;i++) swap(a[t][i],a[r][i]);//在考虑第c列的情况下，将初始行数r的所有数与最大行数t的所有数进行交换
for(int i=n;i>=c;i--) a[r][i]/=a[r][c];//将第r行的第c列的数变成1，该行所有数都要进行初等行列式变换
for(int i=r+1;i<n;i++){ //从第r+1行起到第n行 的所有第c列的系数全部通过初等行列式变换成0
if(fabs(a[i][c])>eps){
for(int j=n;j>=c;j--){
a[i][j]-=a[r][j]*a[i][c];
}
}
}
r++;
}
if(r<n){
for(int i=r;i<n;i++){
if(fabs(a[i][n])>eps) return 2;//无解
}
return 1; //无穷多组解
}
//唯一解
for(int i=n-1;i>=0;i--){ //反向从第n行到第1行求Xn的值
for(int j=i+1;j<n;j++){
a[i][n]-=a[i][j]*a[j][n];
}
}
return 0;
}
int main(){
scanf("%d",&n);
for(int i=0;i<n;i++){
for(int j=0;j<n+1;j++){
scanf("%lf",&a[i][j]);
}
}
int t=guass();
if(t==0) for(int i=0;i<n;i++) printf("%.2lf\n",a[i][n]); //有唯一解并输出唯一解
else if(t==1) puts("Infinite group solutions"); //有无穷多组解
else puts("No solution"); //无解
return 0;
}
如果大家看懂了第一题的做法，可以试着解决一下第二题，是第一题的简单变式。
第二题：题目来源：acwing里面的：高斯消元解异或方程组
#include<bits/stdc++.h>
using namespace std;
const int N=105;
int n;
int a[N][N];
int gauss(){
int r,c;
for(r=c=0;c<n;c++){
int t=r;
for(int i=r;i<n;i++){
if(a[i][c]){
t=i;break;
}
}
if(!a[t][c]) continue;
for(int i=c;i<=n;i++) swap(a[t][i],a[r][i]);
for(int i=r+1;i<n;i++){
if(a[i][c]){
for(int j=c;j<=n;j++){
a[i][j]^=a[t][j];
}
}
}
r++;
}
if(r<n){
for(int i=r;i<n;i++){
if(!a[i][n]) return 2;
}
return 1;
}
for(int i=n-1;i>=0;i--){
for(int j=i+1;j<n;j++){
a[i][n]^=a[i][j]*a[j][n];
}
}
return 0;
}
int main(){
cin>>n;
for(int i=0;i<n;i++){
for(int j=0;j<=n;j++){
cin>>a[i][j];
}
}
int res=gauss();
if(res==0){
for(int i=0;i<n;i++) cout<<a[i][n]<<"\n";
}
else if(res==1) cout<<"Infinite group solutions"<<"\n";
else cout<<"No solution"<<"\n";
}
**现在大家应该都对高斯消元有了一定的理解了。下面直接上两大题给大家练练手，第一时间没有解出不用气馁，这里一道是2021年ICPC济南的A题，一道是2020-2021 ACM-ICPC Asia Seoul Regional Contest ----J题 **
这题说白了就是算每一列的自由元有多少个.假设第i列的自由元有Xi个，则每一列对答案的贡献就是2^(Xi)，最终的答案就是每一列对答案的贡献的乘积，数据有点大，在算每一列对自由元的贡献时记得加上快速幂节省时间。 也可以计算先计算出所有的自由元个数num，答案就是2^(num),记得答案要取模。
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const ll mod=998244353;
const int N=205;
int n;
ll a[N][N],b[N][N],d[N][N];
ll ksm(ll a,ll b){
ll res=1;
while(b){
if(b&1) res=res*a%mod;
a=a*a%mod;
b>>=1;
}
return res%mod;
}
ll gauss(){
int r,c,num=0;
for(r=c=1;c<=n;c++){
int t=r;
for(int i=r+1;i<=n;i++){
if(d[i][c]){
t=i;break;
}
}
if(t!=r) for(int i=c;i<=n;i++) swap(d[t][i],d[r][i]);
if(!d[r][c]){
num++;
continue;
}
for(int i=r+1;i<=n;i++){
if(d[i][c]){
for(int j=c;j<=n;j++){
d[i][j]^=d[r][j];
}
}
}
r++;
}
return num;
}
int main(){
ios::sync_with_stdio(0);
cin>>n;
for(int i=1;i<=n;i++){
for(int j=1;j<=n;j++){
cin>>a[i][j];
}
}
for(int i=1;i<=n;i++){
for(int j=1;j<=n;j++){
cin>>b[i][j];
}
}
ll ans=0;
for(int i=1;i<=n;i++){
for(int j=1;j<=n;j++){
for(int k=1;k<=n;k++){
d[j][k]=a[j][k];
}
d[j][j]^=b[j][i];
}
ans+=gauss();
}
ll res=1;
res=ksm(2,ans)%mod;
cout<<res<<endl;
return 0;
}
2020-2021 ACM-ICPC Asia Seoul Regional Contest J-Switches
{
for(int i=1;i<=n;i++){
for(int j=1;j<=n;j++){
a[i][j]=c[i][j];
}
}
}
int gauss(){
int c,r;
for(r=c=1;c<=n;c++){
int t=r;
for(int i=r;i<=n;i++){
if(a[i][c]){
t=i;break;
}
}
if(!a[t][c]) continue;
for(int i=c;i<=n*2;i++) swap(a[t][i],a[r][i]);
for(int i=1;i<=n;i++){
if(a[i][c]&&i!=r){
for(int j=c;j<=n*2;j++){
a[i][j]^=a[r][j];
}
}
}
r++;
}
if(r<=n) return 1;
else return 0;
}
int main(){
scanf("%d",&n);
for(int j=1;j<=n;j++){
for(int i=1;i<=n;i++){
scanf("%d",&a[i][j]);
c[i][j]=a[i][j];
if(i==j) ch[i][j]=1;
}
a[j][n+j]=1; //灯亮
}
int res=gauss();
if(res==1) puts("-1");
else{
for(int j=n+1;j<=n*2;j++){
for(int i=1;i<=n;i++){
if(a[i][j]) printf("%d ",i);
}
printf("\n");
}
}
return 0;
}
bitset优化写法
#include<bits/stdc++.h>
#include<bitset>
using namespace std;
const int N=505;
bitset<N*2>a[N];
bitset<N>b,inv_a[N],ans;
int n;
bool gauss(){
int c;
for(int c=1;c<=n;c++){
int t=c;
for(int i=c;i<=n;i++){
if(a[i][c]){
t=i;
break;
}
}
if(a[t][c]==0) return false;
swap(a[t],a[c]);
for(int i=1;i<=n;i++){
if(a[i][c]&&i!=c){
a[i]^=a[c];
}
}
}
return true;
}
void mul(bitset<N>& a,bitset<N> b[],bitset<N>& c){
for(int i=1;i<=n;i++){
a[i]=(b[i]&c).count()%2;
}
}
int main(){
scanf("%d",&n);
for(int j=1;j<=n;j++){
for(int i=1;i<=n;i++){
int x;
scanf("%d",&x);
if(x) a[i][j]=1;
}
a[j][j+n]=1;
}
if(gauss()==0) puts("-1");
else{
for(int i=1;i<=n;i++){
for(int j=1;j<=n;j++){
inv_a[i][j]=a[i][j+n];
}
}
for(int i=1;i<=n;i++){
b.reset();
b[i]=1;
ans.reset();
mul(ans,inv_a,b);
for(int i=1;i<=n;i++){
if(ans[i]) printf("%d ",i);
}
printf("\n");
}
}
}
2万+
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
