高斯消元学习笔记
高斯消元学习笔记
——by sunzz3183
介绍
高斯消元是一种求解线性方程组的方法。线性方程组就是 \(m\) 个 \(n\) 元一次方程。如：
求解
线性方程组可以写成一个由系数和常数组成的 \(m\) 行 \(n+1\) 列的增广矩阵。如：
增广矩阵有三个性质：
-
任意交换矩阵的两行或两列，矩阵不变；
-
矩阵一行乘上一个非零数，矩阵不变；
-
把其中一行的若干倍加到另一行上，矩阵不变。
那我们运用矩阵的各种性质，来将矩阵消成对角线上的元素为 \(1\)，并且除了第 \(n+1\) 列其余元素均为 \(0\) 的矩阵，这样我们就很容易的得出每个未知数的值：分别是从上到下第 \(n+1\) 列的值（因为这时候第 \(i\) 行的第 \(i\) 个未知数系数都为 \(1\)）。
那么怎么消呢？
高斯约旦消元法
以上面的矩阵为例子：
明确我们的目的：把矩阵消成对角线为 \(1\)，除了第 \(n+1\) 列其余元素都为 \(0\)。
若有一列全为 \(0\)，那么肯定有第 \(i\) 行第 \(i\) 列消不成 \(1\)，此时无解（少了一个方程，未知数就无法确定了）。
知道了这个，我们就可以对这个矩阵进行初步判定：
for(int i=1,j,l;i<=n;i++){
for(j=i+1,t=i;j<=n;j++)
if(fabs(a[j][i])>fabs(a[t][i]))t=j;//选绝对值最大的作用一会讲
for(j=1;j<=n+1;j++)swap(a[i][j],a[t][j]);
if(fabs(a[i][i])<eps)return puts("No Solution"),0;
}
然后利用性质 \(2\) 让第 \(i\) 行第 \(i\) 列系数化为 \(1\)。
for(j=1,k=a[i][i];j<=n+1;j++)a[i][j]/=k;
最后，利用性质 \(3\) 将其他行的系数化为 \(0\)
for(j=1;j<=n;j++)
for(l=1,k=a[j][i];i!=j&&l<=n+1;l++)
a[j][l]-=k*a[i][l];
然后就完成了！
最后第 \(n+1\) 列就是答案了。
完整代码
#include<bits/stdc++.h>
// #define int long long
#define eps 0.000000001
using namespace std;
inline int read(){
char ch=getchar();int x=0;bool f=1;
while(ch<'0'||'9'<ch){if(ch=='-')f=0;ch=getchar();}
while('0'<=ch&&ch<='9'){x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}
return f?x:-x;
}
const int N=1e2+3;
int n,t=1;
double a[N][N],k;
signed main(){
// freopen(".in","r",stdin);
// freopen(".out","w",stdout);
n=read();
for(int i=1;i<=n;i++)
for(int j=1;j<=n+1;j++)
scanf("%lf",&a[i][j]);
for(int i=1,j,l;i<=n;i++){
for(j=i+1,t=i;j<=n;j++)
if(fabs(a[j][i])>fabs(a[t][i]))t=j;
if(fabs(a[t][i])<eps)return puts("No Solution"),0;
for(j=1;j<=n+1;j++)swap(a[i][j],a[t][j]);
for(j=1,k=a[i][i];j<=n+1;j++)a[i][j]/=k;
for(j=1;j<=n;j++)
for(l=1,k=a[j][i];i!=j&&l<=n+1;l++)
a[j][l]-=k*a[i][l];
}
for(int i=1;i<=n;i++)
printf("%.2f\n",a[i][n+1]);
return 0;
}
遗留问题
\(Q\)：为什么要选绝对值最大的数？
\(A\)：在消元过程中，为了避免出现除数为 \(0\) 的情况，需要在每一行中找到一个绝对值最大的非零数作为主元，然后将该行与其他行进行消元。这样可以避免出现较大的误差，提高计算精度。IEEE 754 浮点数标准中数越接近 \(0\) 越精确。
例题
1.模板
数据较水，正常就能过。
2. 进阶
区别：数据更强，并加入了判断无解和多解。
判断方法：如果存在一行全为 \(0\) 显然就是多解，如果存在一行系数全为 \(0\) 常数不为 \(0\)，那么是无解。
注意：无解的优先级比多解高，多解的优先级比唯一解高。
3. 异或消元
异或就是不进位的加法和减法，正常做就行。
优化：可以用一个数表示二进制数字。
代码：
#include<bits/stdc++.h>
// #define int long long
using namespace std;
inline int read(){
char ch=getchar();int x=0;bool f=1;
while(ch<'0'||'9'<ch){if(ch=='-')f=0;ch=getchar();}
while('0'<=ch&&ch<='9'){x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}
return f?x:-x;
}
const int N=31;
int T,n,t=1,ans;
int a[N];
int yudan(){
for(int i=1,j,k;i<=n;i++){
for(j=i+1,t=i;j<=n;j++)
if(a[j]>a[t])t=j;//找主元最大的，实际上就是找a[j][i]==1
swap(a[i],a[t]);
if(!a[i])return 1<<(n-i+1);
if(a[i]==1)return 0;
k=floor(log2(a[i]));
for(j=1;(a[i]>>k&1)&&j<=n;j++)
if(j!=i&&(a[j]>>k&1))a[j]^=a[i];
}
return 1;
}
signed main(){
// freopen(".in","r",stdin);
// freopen(".out","w",stdout);
T=read();
while(T--){
n=read();
for(int i=1;i<=n;i++)a[i]=read(),a[i]|=1<<i;
for(int i=1;i<=n;i++)a[i]^=read();
int u=read(),v=read();
while(u&&v)
a[v]|=1<<u,u=read(),v=read();
int ans=yudan();
if(!ans)puts("Oh,it's impossible~!!");
else printf("%lld\n",ans);
}
return 0;
}
