c++ - Eigen 与 OpenMP : No parallelisation due to false sharing and thread overhead

Question

系统规范:

1. 英特尔至强 E7-v3 处理器(4 个插槽，16 个内核/插槽，2线程/核心)
2. 使用 Eigen 系列和 C++

以下是代码片段的串行实现:

Eigen::VectorXd get_Row(const int j, const int nColStart, const int nCols) {

    Eigen::VectorXd row(nCols);
    for (int k=0; k<nCols; ++k) {
        row(k) = get_Matrix_Entry(j,k+nColStart);
    }

} 

double get_Matrix_Entry(int x , int y){
    return exp(-(x-y)*(x-y));
} 

我需要并行化 get_Row 部分，因为 nCols 可以大到 10^6，因此，我尝试了某些技术:

1. 朴素并行化:

   Eigen::VectorXd get_Row(const int j, const int nColStart, const int nCols) {  
       Eigen::VectorXd row(nCols);
   
        #pragma omp parallel for schedule(static,8)    
        for (int k=0; k<nCols; ++k) {
             row(k)    =   get_Matrix_Entry(j,k+nColStart);
   
        return row;
   }
   

2. 露天采矿:

   Eigen::VectorXd get_Row(const int j, const int nColStart, const int nCols) { 
       int vec_len = 8;
       Eigen::VectorXd row(nCols) ;
       int i,cols;
       cols=nCols;
       int rem = cols%vec_len;
       if(rem!=0)
           cols-=rem;
   
       #pragma omp parallel for    
       for(int ii=0;ii<cols; ii+=vec_len){
            for(i=ii;i<ii+vec_len;i++){
                row(i) = get_Matrix_Entry(j,i+nColStart);
            }
       }
   
       for(int jj=i; jj<nCols;jj++)
           row(jj) = get_Matrix_Entry(j,jj+nColStart);
   
       return row;
   }
   

3. 来自互联网的某个地方以避免虚假分享:

   Eigen::VectorXd get_Row(const int j, const int nColStart, const int nCols) {
       int cache_line_size=8;
       Eigen::MatrixXd row_m(nCols,cache_line_size);
   
       #pragma omp parallel for schedule(static,1)
       for (int k=0; k<nCols; ++k) 
           row_m(k,0)  =   get_Matrix_Entry(j,k+nColStart);
   
       Eigen::VectorXd row(nCols); 
       row = row_m.block(0,0,nCols,1);
   
      return row;
   
   }
   

输出:

以上技术都没有帮助减少为大型 nCols 执行 get_row 所花费的时间，这意味着 naice 并行化与其他技术的工作方式相似(尽管串行更好)，有什么建议或方法可以帮助缩短时间？

正如用户 Avi Ginsburg 所提到的，我提到了一些其他系统细节:

• g++(GCC) 是 4.4.7 版的编译器
• Eigen 库版本为 3.3.2
• 使用的编译器标志:“-c -fopenmp -Wall -march=native -O3 -funroll-all-loops -ffast-math -ffinite-math-only -I header”，这里的 header 是包含 Eigen 的文件夹。<

• gcc -march=native -Q --help=target->的输出(只提到了一些标志的描述):

  -mavx [启用]

  -mfancy-math-387 [启用]

  -mfma [禁用]

  -march=核心2

有关标志的完整说明，请参阅 this .

Answer 1

尝试将您的函数重写为单个表达式并让 Eigen 向量化自身，即:

Eigen::VectorXd get_Row(const int j, const int nColStart, const int nCols) {

    Eigen::VectorXd row(nCols);

    row = (-( Eigen::VectorXd::LinSpaced(nCols, nColStart, nColStart + nCols - 1).array()
                      - double(j)).square()).exp().matrix();

    return row;
}

确保在编译时使用-mavx 和-mfma(或-march=native)。在 i7 上给我 x4 加速(我知道你在谈论尝试使用 64/128 线程，但这是单线程)。

您可以通过将计算分成多个部分来启用 openmp 以进一步提高速度:

Eigen::VectorXd get_Row_omp(const int j, const int nColStart, const int nCols) {

    Eigen::VectorXd row(nCols);

#pragma omp parallel
    {
        int num_threads = omp_get_num_threads();
        int tid = omp_get_thread_num();
        int n_per_thread = nCols / num_threads;
        if ((n_per_thread * num_threads < nCols)) n_per_thread++;
        int start = tid * n_per_thread;
        int len = n_per_thread;
        if (tid + 1 == num_threads) len = nCols - start;

        if(start < nCols)
            row.segment(start, len) = (-(Eigen::VectorXd::LinSpaced(len,
                               nColStart + start, nColStart + start + len - 1)
                            .array() - double(j)).square()).exp().matrix();

    }
    return row;

}

对于我(4 核)，我在计算 10^8 个元素时获得额外的 ~x3.3 加速，但预计 10^6 和/或 64/128 核的速度会更低(对核数进行归一化，当然)。

编辑

我没有进行任何检查以确保 OMP 线程没有越界并且我在串行版本的 Eigen::VectorXd::LinSpaced 中混淆了第二个和第三个参数。这可能是您遇到的任何错误的原因。此外，我在此处粘贴了用于测试的代码。我用 g++ -std=c++11 -fopenmp -march=native -O3编译，适应你的需要。

#include <Eigen/Core>
#include <iostream>
#include <omp.h>


double get_Matrix_Entry(int x, int y) {
        return exp(-(x - y)*(x - y));
}

Eigen::VectorXd get_RowOld(const int j, const int nColStart, const int nCols) {

        Eigen::VectorXd row(nCols);
        for (int k = 0; k<nCols; ++k) {
                row(k) = get_Matrix_Entry(j, k + nColStart);
        }
        return row;
}


Eigen::VectorXd get_Row(const int j, const int nColStart, const int nCols) {

        Eigen::VectorXd row(nCols);

        row = (-( Eigen::VectorXd::LinSpaced(nCols, nColStart, nColStart + nCols - 1).array() - double(j)).square()).exp().matrix();

        return row;
}

Eigen::VectorXd get_Row_omp(const int j, const int nColStart, const int nCols) {

        Eigen::VectorXd row(nCols);

#pragma omp parallel
        {
                int num_threads = omp_get_num_threads();
                int tid = omp_get_thread_num();
                int n_per_thread = nCols / num_threads;
                if ((n_per_thread * num_threads < nCols)) n_per_thread++;
                int start = tid * n_per_thread;
                int len = n_per_thread;
                if (tid + 1 == num_threads) len = nCols - start;


#pragma omp critical
{
        std::cout << tid << "/" << num_threads << "\t" << n_per_thread << "\t" << start <<
                                                         "\t" << len << "\t" << start+len << "\n\n";
}

                if(start < nCols)
                        row.segment(start, len) = (-(Eigen::VectorXd::LinSpaced(len, nColStart + start, nColStart + start + len - 1).array() - double(j)).square()).exp().matrix();

        }
        return row;
}

int main()
{
        std::cout << EIGEN_WORLD_VERSION << '.' << EIGEN_MAJOR_VERSION << '.' << EIGEN_MINOR_VERSION << '\n';
        volatile int b = 3;
        int sz = 6553600;
        sz = 16;
        b = 6553500;
        b = 3;
        {
                auto beg = omp_get_wtime();
                auto r = get_RowOld(5, b, sz);
                auto end = omp_get_wtime();
                auto diff = end - beg;
                std::cout << r.rows() << "\t" << r.cols() << "\n";
//              std::cout << r.transpose() << "\n";
                std::cout << "Old: " << r.mean() << "\n" << diff << "\n\n";

                beg = omp_get_wtime();
                auto r2 = get_Row(5, b, sz);
                end = omp_get_wtime();
                diff = end - beg;
                std::cout << r2.rows() << "\t" << r2.cols() << "\n";
//              std::cout << r2.transpose() << "\n";
                std::cout << "Eigen:         " << (r2-r).cwiseAbs().sum() << "\t" << (r-r2).cwiseAbs().mean() << "\n" << diff << "\n\n";

                auto omp_beg = omp_get_wtime();
                auto r3 = get_Row_omp(5, b, sz);
                auto omp_end = omp_get_wtime();
                auto omp_diff = omp_end - omp_beg;
                std::cout << r3.rows() << "\t" << r3.cols() << "\n";
//              std::cout << r3.transpose() << "\n";
                std::cout << "OMP and Eigen: " << (r3-r).cwiseAbs().sum() << "\t" << (r - r3).cwiseAbs().mean() << "\n" << omp_diff << "\n";
        }

        return 0;

}