高斯消元 矩阵求秩

杭电zhu and 772002 为例,看高斯消元如何使用

矩阵的秩等于矩阵线性无关的列的数量,自由度 = 矩阵列数 - 矩阵的秩。摘自贴吧。

线性代数,忘了。好像是通过行变换,得到一个上三角矩阵,求出某一变元的解,然后回迭求解。


下面是个例子:

#define _CRT_SECURE_ON_WARNINGS  
#include <iostream>  
#include <cstring>  
#include <algorithm>  
using namespace std;
#define MOD 1000000007  
const int maxn = 2017; //素数最大不超过2000  
typedef int Matrix[maxn][maxn];
int Rank(Matrix A, int m, int n) { //增广矩阵A,m行数(方程式个数),n列数。均从0开始
	int i = 0, j = 0, k, r, u;
	while (i < m && j < n) { //显然j和列有关
		r = i; //选主行  
		for (k = i; k < m; ++k) {
			if (A[k][j]) {
				r = k;
				break;
			}
		}
		if (A[r][j]) {
			if (r != i)
				for (k = 0; k <= n; ++k)
					swap(A[r][k], A[i][k]); //行交换,将第i行换为A[i][j]不为0
			for (u = i + 1; u < m; ++u) //将第j列第i行一下的值全变为0。当然想改变A[u][j],u这一行都要做相同的改变(即全异或A[i][k])
				if (A[u][j])
					for (k = i; k <= n; ++k)
						A[u][k] ^= A[i][k]; //若不是异或,普通的矩阵改为减除
			++i; //移动主行
		}
		++j;//移动列
	}
	return i;//直接返回秩了,没有回带求解
}
Matrix A;
bool vis[maxn];
int prime[304]; //存素数  
int getPrime(int n) { //O(n)线性素数筛  
	memset(vis, 1, sizeof vis);
	int cnt = 0; //从0开始存,若换用注释中的,则从1开始存  
	for (int i = 2; i <= n; ++i) {
		if (vis[i]) {
			prime[cnt++] = i;//p[++cnt] = i;  
		}
		for (int j = 0; j < cnt && (i * prime[j] <= n); ++j) { //j = 1; j <= cnt   
			vis[i * prime[j]] = 0;
			if (i % prime[j] == 0) break;
		}
	}
	return cnt; //返回2:n中素数的个数  
}
long long pow_mod(int n) {
	long long b = 2, res = 1;
	while (n) {
		if (n & 1) res = (res * b) % MOD;
		b = (b * b) % MOD;
		n >>= 1;
	}
	return res;
}
int main()
{
	ios::sync_with_stdio(false);
	int m = getPrime(2000); //303  
	int T, kase = 0;
	cin >> T;
	while (T--) {
		int n, maxp = 0; //maxp 存最大素数的下标  
		long long x;
		cin >> n;
		memset(A, 0, sizeof A); //每一列表示一个数的的质因子分解
		for (int i = 0; i < n; ++i) {
			cin >> x;
			for (int j = 0; j < m; ++j)
				while (x % prime[j] == 0) {
					maxp = max(maxp, j);
					x /= prime[j];
					A[j][i] ^= 1; //第i个元素是否含有奇数个质因子prime[j]  
				}
		} //求出增广矩阵  
		int r = Rank(A, maxp + 1, n); //求秩  
		long long ans = pow_mod(n - r); //自由度 = 矩阵列数 - 矩阵的秩  
		ans = (ans - 1 + MOD) % MOD; //去掉全不选的一种情况  
		cout << "Case #" << ++kase << ":" << endl << ans << endl;
	}
	return 0;
}


POJ1830 套高斯消元的板子,加了一个判断无果的判定,会计算机实现高斯消元后,问题的关键就在于写出方程组,为何要求自由元。

#define _CRT_SECURE_ON_WARNINGS
#include <iostream>
#include <cstring>
#include <algorithm>
using namespace std;
#define MOD 1000000007
const int maxn = 50;
typedef int Matrix[maxn][maxn];
int Rank(Matrix A, int m, int n) {
	int i = 0, j = 0, k, r, u;
	while (i < m && j < n) { //第n列为矩阵右边的值
		r = i; //选主行
		for (k = i; k < m; ++k) {
			if (A[k][j]) { //向下选,定主行
				r = k;
				break;
			}
		}
		if (A[r][j]) {
			if (r != i)//主行上移,i行才是真正的主子
				for (k = 0; k <= n; ++k)
					swap(A[r][k], A[i][k]);
			for (u = i + 1; u < m; ++u)
				if (A[u][j])
					for (k = i; k <= n; ++k)
						A[u][k] ^= A[i][k];
			++i;
		}
		++j;
	}
	for (int k = i; k < m; ++k)
		if (A[k][n] != 0) return -1;//左为0,等号右不为0
	return i;
}
Matrix A;
long long pow_mod(int n) {
	long long b = 2, res = 1;
	while (n) {
		if (n & 1) res = (res * b) % MOD;
		b = (b * b) % MOD;
		n >>= 1;
	}
	return res;
}
int main()
{
	ios::sync_with_stdio(false);
	int T, kase = 0;
	cin >> T;
	while (T--) {
		int n;
		cin >> n;
		memset(A, 0, sizeof A);
		for (int i = 0; i < n; ++i)
			cin >> A[i][n];
		int tmp;
		for (int i = 0; i < n; ++i) {
			cin >> tmp;
			A[i][n] ^= tmp;
			A[i][i] = 1;
		}
		int i, j;
		while (cin >> i >> j) {
			if (i == 0 && j == 0) break;
			A[j - 1][i - 1] = 1; //第i个开关影响第j个
		}
		int r = Rank(A, n, n); //求秩
		if (r == -1) {
			puts("Oh,it's impossible~!!");
			continue;
		}
		long long ans = pow_mod(n - r); //自由度 = 矩阵列数 - 矩阵的秩
		cout << ans << endl;
	}
	return 0;
}






  • 0
    点赞
  • 13
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
以下是一个用C语言实现高斯消元复数逆矩阵的示例代码: ```c #include <stdio.h> #include <stdlib.h> #include <math.h> #define MAX_SIZE 100 void print_matrix(int n, double A[][2*MAX_SIZE]) { int i, j; for (i = 0; i < n; i++) { for (j = 0; j < 2*n; j++) { printf("%lf ", A[i][j]); } printf("\n"); } } void swap_rows(int n, double A[][2*MAX_SIZE], int i, int j) { int k; for (k = 0; k < 2*n; k++) { double temp = A[i][k]; A[i][k] = A[j][k]; A[j][k] = temp; } } void scale_row(int n, double A[][2*MAX_SIZE], int i, double scalar) { int j; for (j = 0; j < 2*n; j++) { A[i][j] *= scalar; } } void add_rows(int n, double A[][2*MAX_SIZE], int i, int j, double scalar) { int k; for (k = 0; k < 2*n; k++) { A[i][k] += scalar * A[j][k]; } } int find_pivot(int n, double A[][2*MAX_SIZE], int k) { int i, max_row = k; double max_val = fabs(A[k][k]); for (i = k+1; i < n; i++) { double val = fabs(A[i][k]); if (val > max_val) { max_val = val; max_row = i; } } return max_row; } void gaussian_elimination(int n, double A[][2*MAX_SIZE]) { int i, j, k; for (k = 0; k < n; k++) { int pivot_row = find_pivot(n, A, k); if (pivot_row != k) { swap_rows(n, A, k, pivot_row); } double pivot = A[k][k]; scale_row(n, A, k, 1.0/pivot); for (i = k+1; i < n; i++) { double factor = -A[i][k]; add_rows(n, A, i, k, factor); } } } void backward_substitution(int n, double A[][2*MAX_SIZE]) { int i, j, k; for (k = n-1; k >= 0; k--) { for (i = 0; i < k; i++) { double factor = -A[i][k]; add_rows(n, A, i, k, factor); } } } void invert_matrix(int n, double A[][2*MAX_SIZE]) { int i, j; for (i = 0; i < n; i++) { for (j = n; j < 2*n; j++) { if (i == j-n) { A[i][j] = 1.0; } else { A[i][j] = 0.0; } } } gaussian_elimination(n, A); backward_substitution(n, A); } int main() { int n, i, j; double A[MAX_SIZE][2*MAX_SIZE]; printf("Enter the size of the matrix: "); scanf("%d", &n); printf("Enter the matrix:\n"); for (i = 0; i < n; i++) { for (j = 0; j < n; j++) { printf("A[%d][%d] = ", i, j); scanf("%lf", &A[i][j]); } } invert_matrix(n, A); printf("The inverse matrix is:\n"); for (i = 0; i < n; i++) { for (j = n; j < 2*n; j++) { printf("%lf ", A[i][j]); } printf("\n"); } return 0; } ``` 该程序使用了高斯消元算法,分为三个部分: 1. `gaussian_elimination`:将矩阵进行高斯消元,将其变为上三角矩阵; 2. `backward_substitution`:将上三角矩阵进行回带解,得到逆矩阵; 3. `invert_matrix`:将原矩阵和单位矩阵合并,进行高斯消元和回带解,得到逆矩阵。 该程序支持任意复数矩阵的逆矩阵解,其中复数使用一个长度为2的数组表示,第一个元素表示实部,第二个元素表示虚部。

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值