具有 mkl 后端的特征库的系数数组操作的性能

Question

我正在将一个包含大量系数数组操作的 Matlab 算法移植到 C++，它看起来像这个示例，但通常要复杂得多：

Eigen::Array<double, Dynamic, 1> tx2(12);
tx2 << 1,2,3,4,5,6;
Eigen::Array<double, Dynamic, 1> tx1(12);
tx1 << 7,8,9,10,11,12;
Eigen::Array<double, Dynamic, 1> x = (tx1 + tx2) / 2;

C++ 代码比 Matlab 慢得多（大约 20%）。因此，在下一步中，我尝试打开 Eigen 的英特尔 MKL 实现，这对性能没有任何作用，就像字面上没有任何改进一样。 MKL 是否有可能不改进系数向量运算？有没有办法测试我是否成功链接 MKL？ Eigen::vector 类是否有更快的替代方案？ 提前致谢！

编辑：我在运行 win7 64 位的 i7-3820 上使用 VS 2013。 更长的例子是：

    Array<double, Dynamic, 1> ts = (k2 / (6 * b.pow(3)) + k / b - b / 2) - (k2 / (6 * a.pow(3)) + k / a - a / 2);
    Array<double, Dynamic, 1> tp1 = -2 * r2*(b - a)/ (rp.pow(2));
    Array<double, Dynamic, 1> tp2 = -2 * r2*rp*log(b / a) / rm2;
    Array<double, Dynamic, 1> tp3 = r2*(b.pow(-1) - a.pow (-1)) / 2;
    Array<double, Dynamic, 1> tp4 = 16 * r2.pow(2)*(r2.pow(2) + 1)*log((2 * rp*b - rm2) / (2 * rp*a - rm2)) / (rp.pow(3)*rm2);
    Array<double, Dynamic, 1> tp5 = 16 * r2.pow(3)*((2 * rp*b - rm2).pow(-1) - (2 * rp*a - rm2).pow(-1)) / rp.pow(3);
    Array<double, Dynamic, 1> tp = tp1 + tp2 + tp3 + tp4 + tp5;
    Array<double, Dynamic, 1> f = (ts + tp) / (2 * ds*ds);

CMakeLists 的相关部分

    set (CMAKE_C_FLAGS "${CMAKE_C_FLAGS} ${OpenMP_C_FLAGS}")
    set (CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} ${OpenMP_CXX_FLAGS}")
    target_link_libraries(MK ${VTK_LIBRARIES} ${Boost_LIBRARIES} mkl_intel_lp64_dll.lib mkl_intel_thread_dll.lib mkl_core_dll.lib libiomp5md.lib)

到目前为止，我只定义了 EIGEN_USE_MKL_ALL。

Answer 1

将对pow(2)、pow(3) 的调用替换为square()、cube()。与pow(-1) 相同，它有利地由一个部门代替。我希望 MatLab 能够为您完成所有这些优化，但在 C++ 中，只有在编译器级别工作才能使此类编译时优化成为可能。

Answer 2

简而言之，如果您有英特尔的 C++ 编译器，请使用它。

我构建了一个 MCVE 来测试这里所做的一些假设。我们要测试

1. MKL 的链接
2. 特征向量化
   1. 加法
   2. 乘法
   3. pow(double)
3. 编译器的效果

使用 Visual Studio 2013。

#include <iostream>

//#define EIGEN_DONT_VECTORIZE

// SSE>2 doesn't affect these tests
#ifndef EIGEN_DONT_VECTORIZE // Not needed with Intel C++ Compiler XE 15.0
    #define EIGEN_VECTORIZE_SSE4_2
    #define EIGEN_VECTORIZE_SSE4_1
    #define EIGEN_VECTORIZE_SSSE3
    #define EIGEN_VECTORIZE_SSE3
#endif

#define EIGEN_USE_MKL_ALL 

#include <Eigen/Core>
#include <ctime>
#include <chrono>

#include <mkl.h>

int main(int argc, char* argv[])
{
    srand(time(NULL));
    std::cout << Eigen::SimdInstructionSetsInUse() << "\n";

    int sz = 32 * 1024 * 1024;

    double dummyAdd, dummyMult, dummyPow;

    // Quick test to show linking worked
    {
        float a[16] = {23.54f};
        float r[16] = {0.f};
        float b = 2.f;

        vsPowx(4, a, b, r);
        std::cout << r[0] << "\n";
    }

    Eigen::ArrayXd v1 = Eigen::ArrayXd::Random(sz);
    Eigen::ArrayXd v2 = Eigen::ArrayXd::Random(sz);
    Eigen::ArrayXd v3 = Eigen::ArrayXd::Random(sz);

    auto startTime = std::chrono::high_resolution_clock::now();
    {
        v3 = v1 + v2;
        dummyAdd = v3.sum();
    }
    auto endTime = std::chrono::high_resolution_clock::now();

    std::cout << "Total Time (addition) " << dummyAdd << " = " <<
        std::chrono::duration_cast<std::chrono::milliseconds>(endTime - startTime).count()
        << " milliseconds.\n";

    startTime = std::chrono::high_resolution_clock::now();
    {
        v1 = v3 * v2;   // 
        dummyMult = v1.sum();
    }
    endTime = std::chrono::high_resolution_clock::now();

    std::cout << "Total Time (multiplication) " << dummyMult << " = " <<
        std::chrono::duration_cast<std::chrono::milliseconds>(endTime - startTime).count()
        << " milliseconds.\n";

    startTime = std::chrono::high_resolution_clock::now();
    {
        v3 = v1.pow(3.5);   // 
        dummyPow = v3.sum();
    }
    endTime = std::chrono::high_resolution_clock::now();

    std::cout << "Total Time (pow(3.5)) " << dummyPow << " = " <<
        std::chrono::duration_cast<std::chrono::milliseconds>(endTime - startTime).count()
        << " milliseconds.\n";

    return 0;

}

然后我使用 cl（VS 的编译器）和 Intel C++ Compiler XE 15.0 进行编译，包括和不包括 EIGEN_DONT_VECTORIZE 和 EIGEN_USE_MKL_ALL。我为这些测试编译了 without omp。我得到了 (i5 3470) 有趣的结果。对于 cl，我认为 MKL 是否链接没有区别，但有轻微的

无
554.132
总时间（加法）-2006.37 = 130 毫秒。
总时间（乘法）1.11832e+007 = 137 毫秒。
总时间 (pow(3.5)) -1.#IND = 1730 毫秒。

和

上交所、上交所2
554.132
总时间（加法）-689.959 = 86 毫秒。
总时间（乘法）1.1175e+007 = 87 毫秒。
总时间 (pow(3.5)) -1.#IND = 1695 毫秒。

所以我们看到加法和乘法似乎是矢量化的，但pow 不受 MKL 影响。

英特尔编译器在行为上表现出类似的结果，但使用pow 时效果更好。

无
554.132
总时间（加法）7594.98 = 96 毫秒。
总时间（乘法）1.11818e+007 = 94 毫秒。
总时间 (pow(3.5)) -1.#IND = 921 毫秒。

和

SSE、SSE2、SSE3、SSSE3、SSE4.1、SSE4.2
554.132
总时间（加法）-1953.37 = 87 毫秒。
总时间（乘法）1.11796e+007 = 87 毫秒。
总时间 (pow(3.5)) -1.#IND = 838 毫秒。

没有 EIGEN_USE_MKL_ALL 和

SSE、SSE2、SSE3、SSSE3、SSE4.1、SSE4.2
554.132
总时间（加法）1512.55 = 87 毫秒。
总时间（乘法）1.11759e+007 = 89 毫秒。
总时间 (pow(3.5)) -1.#IND = 843 毫秒。

使用 EIGEN_USE_MKL_ALL。

我可以理解英特尔的编译器倾向于对代码进行超级优化，直至与 MKL 的性能相匹配。我本来希望看到 cl 性能有所不同。最重要的是，如果您需要更好的性能，请使用英特尔 C++ 编译器。