想法是将两个矩阵相乘。并使用本征进行相同的乘法,然后检查结果是否相同。
在下面的代码中,N = 2
返回same thing
,但N = 1000
返回NOT same thing
。为什么?
#include <cstdlib>
#include <iostream>
#include <Eigen/Dense>
using namespace std;
using namespace Eigen;
const int N = 1000;
void mult_matrix(double x[N][N], double y[N][N], double z[N][N]) {
int rows = N;
int cols = N;
for (int i = 0; i < rows; i++)
for (int j = 0; j < cols; j++)
for (int k = 0; k < cols; k++)
z[i][j] += x[i][k] * y[k][j];
}
void check(double *x, double *y, double *z) {
Matrix<double, Dynamic, Dynamic, RowMajor> m =
Matrix<double, Dynamic, Dynamic, RowMajor>::Map(x, N, N) *
Matrix<double, Dynamic, Dynamic, RowMajor>::Map(y, N, N);
cout << m(0, 0) << endl;
cout << Matrix<double, Dynamic, Dynamic, RowMajor>::Map(z, N, N)(0, 0) << endl;
if (m == Matrix<double, Dynamic, Dynamic, RowMajor>::Map(z, N, N))
cout << "same thing" << endl;
else
cout << "NOT same thing" << endl;
}
int main() {
double *a = (double*)malloc(N*N*sizeof(double));
double *b = (double*)malloc(N*N*sizeof(double));
double *c = (double*)malloc(N*N*sizeof(double));
Matrix<double, Dynamic, Dynamic, RowMajor>::Map(a, N, N).setRandom();
Matrix<double, Dynamic, Dynamic, RowMajor>::Map(b, N, N).setRandom();
Matrix<double, Dynamic, Dynamic, RowMajor>::Map(c, N, N).setZero();
mult_matrix((double (*)[N])a, (double (*)[N])b, (double (*)[N])c);
check(a, b, c);
}
最佳答案
Eigen提供成员函数 isApprox()
,该函数可用于检查两个矩阵在数值精度的范围内是否相等。
在您的代码中,可以通过将==
替换为isApprox()
来简单地实现这种比较,如下所示:
if (m.isApprox(Matrix<double, Dynamic, Dynamic, RowMajor>::Map(z, N, N)))
cout << "same thing" << endl;
else
cout << "NOT same thing" << endl;
可以将所需的精度作为可选的第二个参数传递给
isApprox()
。如评论中所讨论的,可能总是存在这样的比较可能无法可靠进行的情况。但是,使用Eigen的
isApprox()
或isMuchSmallerThan()
等功能比任何简单的手工解决方案都更有效。