最近,我开始使用Eigen(3.3.1版),在OLS回归的核心,即计算矩阵乘积的逆,对Armadillo运行一个基准测试。我注意到,在使用MKL库编译时,Eigen的运行速度比不使用它时慢。我想知道我的编译说明是否错误。我还试图通过直接调用MKL BLAS和LAPACK例程来实现这一操作,并获得了更快的结果,与犰狳一样快。我无法解释如此糟糕的性能,特别是对于浮点类型。
我编写了以下代码来实现这个基准:
#define ARMA_DONT_USE_WRAPPER
#define ARMA_NO_DEBUG
#include <armadillo>
#define EIGEN_NO_DEBUG
#define EIGEN_NO_STATIC_ASSERT
#define EIGEN_USE_MKL_ALL
#include <Eigen/Dense>
template <typename T>
using Matrix = Eigen::Matrix<T, Eigen::Dynamic, Eigen::Dynamic>;
#ifdef USE_FLOAT
using T = float;
#else
using T = double;
#endif
int main()
{
arma::wall_clock timer;
int niter = 1000000;
int n = 1000;
int k = 20;
arma::Mat<T> Xa = arma::cumsum(arma::randn<arma::Mat<T>>(n, k));
Matrix<T> Xe = Matrix<T>::Map(Xa.memptr(), Xa.n_rows, Xa.n_cols);
// Armadillo compiled with MKL
timer.tic();
for (int i = 0; i < niter; ++i) {
arma::Mat<T> iX2a = (Xa.t() * Xa).i();
}
std::cout << "...Elapsed time: " << timer.toc() << "\n";
// Eigen compiled with MKL
timer.tic();
for (int i = 0; i < niter; ++i) {
Matrix<T> iX2e = (Xe.transpose() * Xe).inverse();
}
std::cout << "...Elapsed time: " << timer.toc() << "\n";*/
// Eigen Matrix with MKL routines
timer.tic();
for (int i = 0; i < niter; ++i) {
Matrix<T> iX2e = Matrix<T>::Zero(k, k);
// first stage => computing square matrix trans(X) * X
#ifdef USE_FLOAT
cblas_ssyrk(CblasColMajor, CblasLower, CblasTrans, k, n, 1.0, &Xe(0,0), n, 0.0, &iX2e(0,0), k);
#else
cblas_dsyrk(CblasColMajor, CblasLower, CblasTrans, k, n, 1.0, &Xe(0,0), n, 0.0, &iX2e(0,0), k);
#endif
// getting upper part
for (int i = 0; i < k; ++i)
for (int j = i + 1; j < k; ++j)
iX2e(i, j) = iX2e(j, i);
// second stage => inverting square matrix
// initializing pivots
int* ipiv = new int[k];
// factorizing matrix
#ifdef USE_FLOAT
LAPACKE_sgetrf(LAPACK_COL_MAJOR, k, k, &iX2e(0,0), k, ipiv);
#else
LAPACKE_dgetrf(LAPACK_COL_MAJOR, k, k, &iX2e(0,0), k, ipiv);
#endif
// computing the matrix inverse
#ifdef USE_FLOAT
LAPACKE_sgetri(LAPACK_COL_MAJOR, k, &iX2e(0,0), k, ipiv);
#else
LAPACKE_dgetri(LAPACK_COL_MAJOR, k, &iX2e(0,0), k, ipiv);
#endif
delete[] ipiv;
}
std::cout << "...Elapsed time: " << timer.toc() << "\n";
}
我编译了这个名为test的文件。cpp与:
我得到以下结果(在Intel(R)Core(TM)i5-3210M CPU上,2.50GHz)
具有MKL的犰狳=>64.0s
具有MKL的特征值=>72.2秒
仅本征值=>68.7秒
纯MKL=>64.9秒
具有MKL的犰狳=>38.2s
具有MKL的特征值=>61.1秒
仅本征值=>42.6s
纯MKL=>38.9秒
注意:我为一个不会使用非常大的矩阵的项目运行这个测试,我不需要在这个级别上进行并行化,我最大的矩阵可能是25列2000行,而且我需要在更高的级别上并行,因此我希望避免任何可能会降低代码速度的嵌套并行。