在 ODE 函数中使用 ODEINT 的 OpenMP答案

【问题标题】：OpenMP with ODEINT in ODE function在 ODE 函数中使用 ODEINT 的 OpenMP
【发布时间】：2020-06-05 04:10:41
【问题描述】：

我正在尝试在内部并行化由 ODEINT 集成的 ODE 函数。

我做了以下小例子

#include <iostream>
#include <chrono>
#include <Eigen/Dense>
#include <omp.h>
#include <boost/numeric/odeint.hpp>
#include <boost/numeric/odeint/external/openmp/openmp.hpp>
#include <boost/numeric/odeint/external/eigen/eigen.hpp>

using namespace boost::numeric::odeint;

class System {
private:
    Eigen::VectorXd _input_data;
public:
    System( Eigen::VectorXd &input_data ) { _input_data = input_data; };
    void operator() ( const Eigen::VectorXd &x , Eigen::VectorXd &dxdt , const double t ) {
        double _sum = 0.;
        #pragma omp parallel for reduction(+:_sum)
        for(int k = 0; k < _input_data.size(); ++k) {
            _sum += _input_data(k);
        };
        dxdt(0) = _sum;
    };
};

int main() {
    omp_set_num_threads(1);
    Eigen::VectorXd input_data = Eigen::VectorXd::Zero(100);
    System ode(input_data);
    runge_kutta_dopri5<Eigen::VectorXd> rk5_stepper;
    Eigen::VectorXd x = Eigen::VectorXd::Zero(1);
    auto start = std::chrono::high_resolution_clock::now();
    size_t steps = integrate_const(rk5_stepper, ode, x, 0., 1., 0.01);
    auto stop = std::chrono::high_resolution_clock::now();
    auto duration = std::chrono::duration_cast<std::chrono::microseconds>(stop - start);
    std::cout << "Execution time: " << duration.count() / 1000000. << " sec" << std::endl;
    return 0;
}

使用 CMakeLists.txt 文件

cmake_minimum_required(VERSION 3.13)

set(CMAKE_C_COMPILER /usr/local/bin/gcc-9)
set(CMAKE_CXX_COMPILER /usr/local/bin/g++-9)

project(ODEINT_OPENMP_TEST)

set(CMAKE_CXX_STANDARD 14)

include_directories("/usr/local/include")

find_package(OpenMP REQUIRED)

add_executable(ODEINT_OPENMP_TEST main.cpp)

target_link_libraries(ODEINT_OPENMP_TEST PRIVATE OpenMP::OpenMP_CXX)

当我尝试通过omp_set_num_threads(N) 使用更多线程时，与仅使用单个线程omp_set_num_threads(1) 相比，该程序始终变慢。选择 N=2，程序会变慢大约 x3（在我的机器上）。

直观地说，右边的函数应该并行运行得更快？我做错了吗？

【问题讨论】：

标签： c++ openmp odeint

【解决方案1】：

首先，您的循环太小，无法从在这个特定示例中使用多线程中受益（在我的机器上按顺序大约需要 3.5 毫秒，在 6 核上使用 6 个线程时需要 1.8 毫秒）。

此外，您的基准测试太短了，您可能会测量到意想不到的影响（缓存、页面错误、处理器频率扩展问题等）。考虑将其置于循环中以减轻大多数影响（如果这在现实条件下有意义的话）。

此外，一些 OpenMP 运行时会在执行并行部分时创建线程。这个操作很慢。由于指令#pragma omp parallel 包含在计时中，您还可以测量线程创建。

这是在我的 6 核机器上的结果，它的大小是大 1000 倍：

1 thread:  2.330 sec
2 threads: 1.212 sec
3 threads: 0.813 sec
4 threads: 0.649 sec
5 threads: 0.532 sec
6 threads: 0.459 sec

6线程加速5.1，不错。

请注意，由于您的循环似乎受内存限制（内存吞吐量不会随着使用的内核数量而扩展），因此您的扩展可能会更差。

【讨论】：