闲聊高斯消元入门
来源：互联网 发布：txt本地听书软件 编辑：程序博客网 时间：2024/06/08 01:22
前言：前一段时间看了看高斯消元，一直想把学习过程写下来，但有点懒=-=现在补上
一、关于高斯消元的题目常见的类型有：
1 . 求唯一整数解
2 . 求唯一浮点数解
3 . 求解得个数，等同于求自由元个数，矩阵的秩。
然后其中一般都会让你判断无解的情况
二、高斯消元的步骤
第一步：将给定矩阵化成阶梯矩阵。
第二步：根据题目的要求 ：
1 . 判断是否无解
2 . 求出解的个数是否唯一
3 . 求出唯一解。
这是整体的思路，细节我就不详细讲了，难度都不大，而且网上都有好文章=-=
三、例题
注：以下代码可能差异很大，因为我也是一边学一边改进的。。
#include<cstdio>#include<cstring>#include<algorithm>using namespace std;const int mod = 7;const int MAXN = 355;/****记得两个数组最好初始化为0****/int mar[MAXN][MAXN];int ans[MAXN];inline int q_pow(int x,int y){ int res = 1; while(y){ if(y&1){ res = res * x % mod; }y>>=1; x = x * x % mod; }return res ;}int gcd(int a,int b){ return b?gcd(b,a%b):a;}inline int lcm(int a,int b){ return a/gcd(a,b)*b;}void debug(int n,int m){ for(int i=0;i<n;i++){ for(int j=0;j<=m;j++){ printf("%d%c",mar[i][j],j!=m?' ':'\n'); } }printf("\n");}/****** n个方程,m个变量******/int gauss(int n,int m){ //debug(n,m);/*****1.转换为阶梯式*******/ for(int col=0,row=0;col<=m&&row<n;col++,row++){ int mx = row; for(int i=row+1;i<n;i++){ if(mar[mx][col]<mar[i][col]){ mx = i; } } if(mar[mx][col]==0){ row--; continue; } for(int i=col;i<=m;i++){ swap(mar[mx][i],mar[row][i]); } //debug(n,m); for(int i=row+1;i<n;i++){ if(mar[i][col]){//debug(n,m); int d = lcm(mar[row][col],mar[i][col]); //printf("d=%d\n",d); int a = d/mar[row][col],b = d/mar[i][col]; for(int j=col;j<=m;j++){//printf("j=%d\n",j); mar[i][j] =( ( mar[row][j]*a - mar[i][j]*b ) % mod + mod ) % mod; } } } } //debug(n,m); //printf("in\n");/******检查是否无解********/int k = 0;/***k为无效的方程****/for(int i=0;i<n;i++){ for(int j=0;j<=m;j++){ if(j!=m&&mar[i][j]){ break; }else if(j==m){ if(i<n&&mar[i][m]){ return 0 ; }else if(!mar[i][m]){ k++; } } }}if(n-k<m){ return 2;}/******从最后一层阶梯向上依次求出未知量的值******/for(int i=m-1;i>=0;i--){ int val = mar[i][m];/***解为浮点数修改处***/ for(int j=m-1;j>i;j--){ val = ( ( val - ans[j] * mar[i][j] ) % mod +mod ) % mod; }ans[i] = (val * q_pow(mar[i][i],mod-2) % mod + mod)%mod;}return 1;}inline int get(char*s){ switch(s[0]){ case 'M':return 1; case 'T':return s[1]=='U'?2:4; case 'W':return 3; case 'F':return 5; case 'S':return s[1]=='A'?6:7; }}int main(){ int n,m; while(scanf("%d%d",&n,&m)==2&&n+m){ memset(mar,0,sizeof(mar)); memset(ans,0,sizeof(ans)); char s1[5],s2[5]; for(int i=0,t,a;i<m;i++){ scanf("%d%s%s",&t,s1,s2); int t1 = get(s1),t2 = get(s2); mar[i][n]=(t2-t1+1+mod)%mod; while(t--){ scanf("%d",&a); a--; mar[i][a]++; mar[i][a]%=mod; } }int k = gauss(m,n); if(k==1){ for(int i=0;i<n;i++){ printf("%d%c",ans[i]<3?ans[i]+7:ans[i],i!=n-1?' ':'\n'); } }else if(k==0){ printf("Inconsistent data.\n"); }else { printf("Multiple solutions.\n"); } }return 0;}类型：求唯一浮点数解
代码：
#include<queue>#include<cstdio>#include<cstring>#include<algorithm>#include<cmath>using namespace std;const int MAXN = 300;double mat[MAXN][MAXN];double ans[MAXN];int gauss(int n,int m){ for(int col=0,row=0;col<=m&&row<n;col++,row++){ int mx = row; for(int i=row+1;i<n;i++){ if(fabs(mat[mx][col])<fabs(mat[i][col])){ mx = i; } } if(mat[mx][col]==0){ row--; continue; }int sign = mat[mx][col]<0?-1:1; for(int i=col;i<=m;i++){ swap(mat[mx][i],mat[row][i]); mat[row][i]*=sign; } for(int i=row+1;i<n;i++){ if(mat[i][col] != 0){ double t = mat[i][col]/mat[row][col]; for(int j=col;j<=m;j++){ mat[i][j] -= mat[col][j]*t ; } } } }for(int i=m-1;i>=0;i--){ double val = 1.0 * mat[i][m]; for(int j=m-1;j>i;j--){ val = ( val - ans[j] * mat[i][j] ); }ans[i] = val / mat[i][i] ;}return 1;}const int X[4]={-1,0,0,1};const int Y[4]={0,1,-1,0};int n,m;char mp[20][20];int id[20][20];;bool is[20][20];bool vis[20][20];queue< pair<int,int> > Q;#define m_p(x,y) make_pair(x,y)#define fi first#define se secondinline int init(){ bool res = false ; while(!Q.empty()){ int sz = Q.size(); while(sz--){ int x = Q.front().fi; int y = Q.front().se; if(mp[x][y]=='$'){ is[x][y]=true; res = true; }Q.pop(); for(int i=0,a,b;i<4;i++){ a = x + X[i]; b = y + Y[i]; if(a>0&&b>0&&a<=n&&b<=m&&mp[a][b]!='#'){ if(!vis[a][b]){ vis[a][b] = true; Q.push(m_p(a,b)); } } } } }int all = 0; for(int i=1;i<=n;i++){ for(int j=1;j<=m;j++){ if(vis[i][j]&&!is[i][j]){ id[i][j]=all++; } } } for(int i=1;i<=n;i++){ for(int j=1;j<=m;j++){ if(vis[i][j]&&!is[i][j]){ int e = id[i][j]; double k = 0.; for(int l=0;l<4;l++){ int a = i+X[l]; int b = j+Y[l]; if(a>0&&b>0&&a<=n&&b<=m&&vis[a][b]){ if(is[a][b]){ k+=1.0; continue; } int& t = id[a][b]; mat[e][t]=-1.0; k++; } } mat[e][e] = k; mat[e][all]= k; } } } if(!res){ return false ; } gauss(all,all); return true;}int main(){ while(scanf("%d%d",&n,&m)==2){ memset(vis,0,sizeof(vis)); memset(is,0,sizeof(is)); memset(mat,0,sizeof(mat)); memset(ans,0,sizeof(ans)); while(!Q.empty())Q.pop(); int x0,y0; for(int i=1;i<=n;i++){ scanf("%s",mp[i]+1); for(int j=1;j<=m;j++){ if(mp[i][j]=='@'){ x0 = i; y0 = j; vis[i][j]=true; Q.push(m_p(i,j)); } } } if(init()){ printf("%.6f\n",ans[id[x0][y0]]); }else { printf("-1\n"); } }return 0;}三、
/***注意： 1.两个数组记得初始化 2.n个方程，m个变量，mar[][0~m-1]存储对应的系数，mar[][m]存储方程的值 3.关于返回值 0：表示无解，阶梯式中出现0,0,0~x(x不为0)的情况。 -1：表示唯一解。 正数：阶梯式中出现的有效方程数目少于变量数目，表示有多组解，（变量数-有效方程数）个无法确定的变元。 4.如果有mod，不会有负数的情况，取逆元进行运算。 5.如果全程浮点数运算，看另一份代码。***/#include<cstdio>#include<cstring>#include<algorithm>using namespace std;const int mod = 2;const int MAXN = 355;int mar[MAXN][MAXN];int ans[MAXN];inline int q_pow(int x,int y){ int res = 1; while(y){ if(y&1){ res = res * x % mod; }y>>=1; x = x * x % mod; }return res ;}int gcd(int a,int b){ return b?gcd(b,a%b):a;}inline int lcm(int a,int b){ return a/gcd(a,b)*b;}void debug(int n,int m){ printf("------------checking------------\n"); for(int i=0;i<n;i++){ for(int j=0;j<=m;j++){ printf("%d%c",mar[i][j],j!=m?' ':'\n'); } }printf("------------checking------------\n");}/****** n个方程,m个变量******/int gauss(int n,int m){/*****1.转换为阶梯式*******/ for(int col=0,row=0;col<=m&&row<n;col++,row++){ int mx = row; for(int i=row+1;i<n;i++){ if(mar[mx][col]<mar[i][col]){ mx = i; } } if(mar[mx][col]==0){ row--; continue; } for(int i=col;i<=m;i++){ swap(mar[mx][i],mar[row][i]); } for(int i=row+1;i<n;i++){ if(mar[i][col]){ int d = lcm(mar[row][col],mar[i][col]); int a = d/mar[row][col],b = d/mar[i][col]; for(int j=col;j<=m;j++){ mar[i][j] =( ( mar[row][j]*a - mar[i][j]*b ) % mod + mod ) % mod;/***mod修改处***/ } } } }/******检查是否无解或者多解********/int k = 0;/***k为无效的方程数目，无效方程为全为0的方程****/for(int i=0;i<n;i++){ for(int j=0;j<=m;j++){ if(j!=m&&mar[i][j]){ break; }else if(j==m){ if(i<n&&mar[i][m]){ return 0 ; }else if(!mar[i][m]){ k++; } } }}/******若有效的方程数目小于变量数，则有（变量数-有效方程数）个无法确定的变元****/if(n-k<m){ return m-n+k;}/******从最后一层阶梯向上依次求出未知量的值******/for(int i=m-1;i>=0;i--){ int val = mar[i][m];/***解为浮点数修改处***/ for(int j=m-1;j>i;j--){ val = ( ( val - ans[j] * mar[i][j] ) % mod +mod ) % mod;/***mod修改处***/ }ans[i] = (val * q_pow(mar[i][i],mod-2) % mod + mod)%mod;/***mod修改处***/}return -1;}int A[MAXN][MAXN];int main(){ int T,t=0;scanf("%d",&T); while(T--){ int n,m; memset(mar,0,sizeof(mar)); memset(ans,0,sizeof(ans)); memset(A,0,sizeof(A)); scanf("%d%d",&n,&m); for(int i=0,t,x;i<m;i++){ scanf("%d",&t); while(t--){ scanf("%d",&x); x--; A[x][i]++; } } printf("Case %d:\n",++t); int q;scanf("%d",&q); while(q--){ for(int i=0;i<n;i++){ for(int j=0;j<m;j++){ mar[i][j]=A[i][j]; } }memset(ans,0,sizeof(ans)); for(int i=0;i<n;i++){ scanf("%d",&mar[i][m]); }//printf("begin------\n"); //debug(n,m); int k = gauss(n,m); if(!k){ printf("0\n"); }else if(k==-1){ printf("1\n"); }else { printf("%I64d\n",1LL<<k); } } }return 0;}四、HDU5833 Zhu and 772002
类型 ： 求解的个数。
代码：
#include<cstdio>#include<cstring>#include<algorithm>#include<cmath>#include<vector>#include<map>#include<set>using namespace std;/******!!!!!!!! Remember !!!; use init() function/***** o(n) ******/namespace prime_table{const int MAX_N = 2000;const long long mod = 2;const long long p = 1e9+7;int all=0;int pr[MAX_N/10+1000];bool isp[MAX_N+10];inline void init(){ all = 0; memset(isp,0,sizeof(isp)); for(int i=2;i<=MAX_N;i++){ if(!isp[i]){ pr[all++] = i; }for(int j=0;j<all;j++){ long long t = 1LL*pr[j]*i ; if(t<=MAX_N){ isp[t] = true; if(i%pr[j]==0)break; }else { break; } } }//printf("all=%d\n",all);return ;}}using namespace prime_table;/**** ×éºÏ*****/long long qpow(long long x,long long y){ long long res =1; while(y){ if(y&1){ res = res * x % p; }x = x * x % p; y >>= 1; }return res ;}long long ni[305];long long pre[305];void iniit(){ ni[0]=0; pre[1] = pre[0] = ni[1]=1; for(int i=2;i<=300;i++){ pre[i] = pre[i-1] * i % p; ni[i] = qpow(i,p-2); }}/***注意： 1.两个数组记得初始化 2.n个方程，m个变量，mar[][0~m-1]存储对应的系数，mar[][m]存储方程的值 3.关于返回值 0：表示无解，阶梯式中出现0,0,0~x(x不为0)的情况。 -1：表示唯一解。 正数：阶梯式中出现的有效方程数目少于变量数目，表示有多组解，（变量数-有效方程数）个无法确定的变元。 4.如果有mod，一般不会有负数的情况，然后取逆元。 5.如果全程浮点数运算，看另一份代码。***/const int MAXN = 355;bool vis[MAXN];long long e[MAXN][MAXN];long long mar[MAXN][MAXN];long long ans[MAXN];inline void oper(int id,long long x){ id--; for(int i=0;i<all&&pr[i]<=x;i++){ int k = 0; while(x%pr[i]==0){ vis[i] = true; e[i][id]=!e[i][id]; x/=pr[i]; } } return ;}inline int work(int n){ int cnt = 0; for(int i=0;i<303;i++){ if(vis[i]){ for(int j=0;j<n;j++){ mar[cnt][j] = e[i][j]; }mar[cnt++][n]=0; } }return cnt;}int gcd(int a,int b){ return b?gcd(b,a%b):a;}inline int lcm(int a,int b){ return a/gcd(a,b)*b;}void debug(int n,int m){ printf("------------checking------------\n"); for(int i=0;i<n;i++,printf("\n")){ for(int j=0;j<=m;j++){ printf("%d%c",mar[i][j],j!=m?' ':'\n'); } }printf("------------checking------------\n");}/****** n个方程,m个变量******/int gauss(int n,int m){//debug(n,m);/*****1.转换为阶梯式*******/ for(int col=0,row=0;col<=m&&row<n;col++,row++){ int mx = row; for(int i=row+1;i<n;i++){ if(mar[mx][col]<mar[i][col]){ mx = i; } } if(mar[mx][col]==0){ row--; continue; } for(int i=col;i<=m;i++){ swap(mar[mx][i],mar[row][i]); } for(int i=row+1;i<n;i++){ if(mar[i][col]){ for(int j=col;j<=m;j++){ mar[i][j] =mar[row][j] ^ mar[i][j] ;/***mod修改处***/ } } } }//debug(n,m);/******检查是否无解或者多解********/int k = 0;/***k为无效的方程数目，无效方程为全为0的方程****/for(int i=0;i<n;i++){ for(int j=0;j<=m;j++){ if(j!=m&&mar[i][j]){ break; }else if(j==m){ if(i<n&&mar[i][m]){ return -1 ; }else if(!mar[i][m]){ k++; } } }}/******若有效的方程数目小于变量数，则有（变量数-有效方程数）个无法确定的变元****/ return m-n+k;}int main(){ //freopen("data.in","r",stdin); //freopen("my_ans.out","w",stdout); init();iniit(); int T,t=0; scanf("%d",&T); while(T--){ int n; scanf("%d",&n); long long x ; memset(vis,0,sizeof(vis)); memset(e,0,sizeof(e)); for(int i=1;i<=n;i++){ scanf("%I64d",&x); oper(i,x); } int res = gauss(work(n),n); //printf("ans = %d\n",res); if(res<=0){ res = 0; }else { res = (int)(( qpow(2,1LL*res)-1 + p)%p); } printf("Case #%d:\n%d\n",++t,res); }return 0;}/***151345 947 1979 1032 1286 1262 879 510 1671 1291 1552 802 1470 1828 15240***/
0 0
- 闲聊高斯消元入门
- 闲聊欧拉函数入门
- 闲聊
- 闲聊
- 闲聊
- 闲聊
- 闲聊
- 闲聊
- 闲聊
- 闲聊
- 闲聊
- 闲聊
- 闲聊
- 闲聊
- 闲聊，与技术无关，我未入门的第一天
- 没事闲聊！
- 闲聊软件工程
- 闲聊软件设计
- java中的字符流
- 面试题
- 解决android SDK 找不到 tools\android.bat
- VC++文件操作之最全篇
- Intersection of Two Arrays
- 闲聊高斯消元入门
- visual studio code的使用
- Android系列--DOM、SAX、Pull解析XML
- Android 开发环境搭建 （Android Studio）
- [2016/8/21][Linux]操作文件和目录
- linux下如何用GDB调试c++程序
- WebLogic中间件
- NGUI批量设置UISprite的Atlas
- QListWidget添加自定义QWidget
