高斯消元同余线性方程的模板

2019-04-14 16:51发布

博主通过POJ 2947 Widget Factory题解总结出来的求解同余线性方程,驾驭扩展矩阵。 题解传送门:http://blog.csdn.net/xf_zhen/article/details/52039622 模板: const int maxn = 305;//数组大小 const int mod = 7; //同模的模 int matrix[305][305],n,m,X[305];//matrix是扩展矩阵 ,n,m为矩阵宽长;X为存储x的值 //bool free_x[305]; int gcd(int a,int b) { return b?gcd(b,a%b):a; } int LCM(int a,int b) { return a*b/gcd(a,b); } /* void Debug(void) { puts(""); int i, j; for (i = 0; i < m; i++) { for (j = 0; j < n + 1; j++) { cout << matrix[i][j] << " "; } cout << endl; } cout << endl; } */ int Guass()//无解返回-1,无数解返回n-k,有解救直接输出并返回0 { int i,j,k,col; memset(X,0,sizeof(X)); //clearall(free_x,1);//把解集清空,所有变量都标为自由变量 //Debug(); for (k = 0,col = 0; k < m && col < n; ++k, ++col) //枚举行列 { int max_r = k;//找到该col列元素绝对值最大的那行与第k行交换.(为了在除法时减小误差) for (i = k + 1; i < m; ++i) { if (abs(matrix[i][col]) > abs(matrix[max_r][col])) max_r = i; } if (max_r != k) //交换 { for (i = col; i < n + 1; ++i) swap(matrix[k][i],matrix[max_r][i]); } //Debug(); if (matrix[k][col] == 0) //如果对应该列都为0,枚举该行的下一列 { k--; continue; } for (i = k + 1; i < m; ++i) //将k后边的col进行初等变换成行阶梯矩阵 { if (matrix[i][col] != 0) { int lcm = LCM(matrix[k][col],matrix[i][col]); int ta = lcm/abs(matrix[i][col]); int tb = lcm/abs(matrix[k][col]); if (matrix[i][col]*matrix[k][col] < 0) tb = -tb; for (j = col; j < n + 1; ++j) { matrix[i][j] =( (ta*matrix[i][j] - tb*matrix[k][j])%mod+mod)%mod; } } } } //Debug();printf("k = %d col = %d matrix = %d ",k,col,matrix[k][col]); // 1. 无解的情况: 化简的增广阵中存在(0, 0, ..., a)这样的行(a != 0). 即R(A) != R(A')无解 for (i = k; i < m; ++i) { if (matrix[i][col] != 0) return -1; } // 2. 无穷解的情况: 在n * (n + 1)的增广阵中出现(0, 0, ..., 0)这样的行,即说明没有形成严格的上三角阵. // 且出现的行数即为自由变元的个数. 即R(A) = R(A') < n //printf("%d %d ",k,n); if (k < n) { //注释处为求多解的自由变量 /*// 首先,自由变元有n - k个,即不确定的变元至少有n - k个. int num = 0,freeidx; for (i = k - 1; i >= 0; --i) { num = 0;// 用于判断该行中的不确定的变元的个数,如果超过1个,则无法求解,它们仍然为不确定的变元. int tmp = matrix[i][n]; // 第i行一定不会是(0, 0, ..., 0)的情况,因为这样的行是在第k行到第m行. // 同样,第i行一定不会是(0, 0, ..., a), a != 0的情况,这样的无解的. for (j = 0; j < n; ++j) { if (matrix[i][j] != 0 && free_x[j]) { num++; freeidx = j; } } if (num > 1) continue; // 无法求解出确定的变元. // 说明就只有一个不确定的变元free_index,那么可以求解出该变元,且该变元是确定的. tmp = matrix[i][n]; for (j = 0; j < n; ++j) { if (matrix[i][j] && j != freeidx) tmp -= matrix[i][j]*X[j]; } X[freeidx] = tmp/matrix[i][freeidx]; free_x[freeidx] = 0; }*/ return n - k; } // 3. 唯一解的情况: 在n * (n + 1)的增广阵中形成严格的上三角阵. // 计算出Xn-1, Xn-2 ... X0. for (i = k - 1; i >= 0; --i) { int tmp = matrix[i][n]; for (j = i + 1; j < n; ++j) { tmp =((tmp- matrix[i][j]*X[j])%mod+mod)%mod; } while(tmp%matrix[i][i]!=0) tmp+=mod; X[i] = (tmp/matrix[i][i])%mod; if(X[i]<3)X[i]+=mod; } return 0; }