我有一个维度为 nxn 的二进制矩阵(零和一)D[][]
,其中 n 很大(大约在 1500 - 2000 左右)。我想在 C 中找到这个矩阵的逆矩阵。
由于我是 C 新手,我从 3 x 3 矩阵开始,然后将其推广到 N x N。这适用于 int 值,但是因为我正在使用二进制1
's 和0
's。在这个实现中,我需要unsigned int
值。
我可以找到许多int
值的解决方案,但我没有遇到任何unsigned int
. 我想在不使用任何外部库(如 blas/lapack)的情况下找到 N x N 二进制矩阵的逆矩阵。如果有人能提供 M x N 矩阵的领先优势,那就太好了。
请注意,我需要矩阵的逆,而不是伪逆。
/* To find the inverse of a matrix using LU decomposition */
/* standard Headers */
#include<math.h>
#include<stdio.h>
int main() {
/* Variable declarations */
int i,j;
unsigned int n,m;
unsigned int rows,cols;
unsigned int D[3][3], d[3], C[3][3];
unsigned int x, s[3][3];
unsigned int y[3];
void LU();
n = 2;
rows=3;cols=3;
/* the matrix to be inverted */
D[0][0] = 1;
D[0][1] = 1;
D[0][2] = 0;
D[1][0] = 0;
D[1][1] = 1;
D[1][2] = 0;
D[2][0] = 1;
D[2][1] = 1;
D[2][2] = 1;
/* Store the matrix value for camparison later.
this is just to check the results, we don't need this
array for the program to work */
for (m = 0; m <= rows-1; m++) {
for (j = 0; j <= cols-1; j++) {
C[m][j] = D[m][j];
}
}
/* Call a sub-function to calculate the LU decomposed matrix. Note that
we pass the two dimensional array [D] to the function and get it back */
LU(D, n);
printf(" \n");
printf("The matrix LU decomposed \n");
for (m = 0; m <= rows-1; m++) {
for (j = 0; j <= cols-1; j++){
printf(" %d \t", D[m][j]);
}
printf("\n");
}
/* TO FIND THE INVERSE */
/* to find the inverse we solve [D][y]=[d] with only one element in
the [d] array put equal to one at a time */
for (m = 0; m <= rows-1; m++) {
d[0] = 0;
d[1] = 0;
d[2] = 0;
d[m] = 1;
for (i = 0; i <= n; i++) {
x = 0;
for (j = 0; j <= i - 1; j++){
x = x + D[i][j] * y[j];
}
y[i] = (d[i] - x);
}
for (i = n; i >= 0; i--) {
x = 0;
for (j = i + 1; j <= n; j++) {
x = x + D[i][j] * s[j][m];
}
s[i][m] = (y[i] - x) / D[i][i];
}
}
/* Print the inverse matrix */
printf("The Inverse Matrix\n");
for (m = 0; m <= rows-1; m++) {
for (j = 0; j <= cols-1; j++){
printf(" %d \t", s[m][j]);
}
printf("\n");
}
/* check that the product of the matrix with its iverse results
is indeed a unit matrix */
printf("The product\n");
for (m = 0; m <= rows-1; m++) {
for (j = 0; j <= cols-1; j++){
x = 0;
for (i = 0; i <= 2; i++) {
x = x + C[m][i] * s[i][j];
}
//printf(" %d %d %f \n", m, j, x);
printf("%d \t",x);
}
printf("\n");
}
return 0;
}
/* The function that calcualtes the LU deomposed matrix.
Note that it receives the matrix as a two dimensional array
of pointers. Any change made to [D] here will also change its
value in the main function. So there is no need of an explicit
"return" statement and the function is of type "void". */
void LU(int (*D)[3][3], int n) {
int i, j, k;
int x;
printf("The matrix \n");
for (j = 0; j <= 2; j++) {
printf(" %d %d %d \n", (*D)[j][0], (*D)[j][1], (*D)[j][2]);
}
for (k = 0; k <= n - 1; k++) {
for (j = k + 1; j <= n; j++) {
x = (*D)[j][k] / (*D)[k][k];
for (i = k; i <= n; i++) {
(*D)[j][i] = (*D)[j][i] - x * (*D)[k][i];
}
(*D)[j][k] = x;
}
}
}
这只是我尝试的一个示例示例,我在逆矩阵中有 -1 值,这是我主要关心的问题。我有 1000 x 1000 的二进制值矩阵,逆矩阵也应该是二进制的。
矩阵:
1 1 0
0 1 0
1 1 1
矩阵 LU 分解:
1 1 0
0 1 0
1 0 1
逆矩阵:
1 -1 0
0 1 0
-1 0 1
产品:
1 0 0
0 1 0
0 0 1