想法是将两个矩阵相乘。并使用本征进行相同的乘法,然后检查结果是否相同。

在下面的代码中,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()等功能比任何简单的手工解决方案都更有效。

10-08 18:31