测量 OpenMP Fork/Join 延迟

Question

由于 MPI-3 具有共享内存并行性的功能，而且它似乎非常适合我的应用程序，我正在认真考虑将我的混合 OpemMP-MPI 代码重写为纯 MPI 实现。

为了将最后一颗钉子钉进棺材，我决定运行一个小程序来测试OpenMPfork/join机制的延迟。这是代码（为英特尔编译器编写）：

void action1(std::vector<double>& t1, std::vector<double>& t2)
{
    #pragma omp parallel for schedule(static) num_threads(std::thread::hardware_concurrency())
    for (auto index = std::size_t{}; index < t1.size(); ++index)
    {
        t1.data()[index] = std::sin(t2.data()[index]) * std::cos(t2.data()[index]);
    }
}

void action2(std::vector<double>& t1, std::vector<double>& t2)
{
    #pragma omp parallel for schedule(static) num_threads(std::thread::hardware_concurrency())
    for (auto index = std::size_t{}; index < t1.size(); ++index)
    {
        t1.data()[index] = t2.data()[index] * std::sin(t2.data()[index]);
    }
}

void action3(std::vector<double>& t1, std::vector<double>& t2)
{
    #pragma omp parallel for schedule(static) num_threads(std::thread::hardware_concurrency())
    for (auto index = std::size_t{}; index < t1.size(); ++index)
    {
        t1.data()[index] = t2.data()[index] * t2.data()[index];
    }
}

void action4(std::vector<double>& t1, std::vector<double>& t2)
{
    #pragma omp parallel for schedule(static) num_threads(std::thread::hardware_concurrency())
    for (auto index = std::size_t{}; index < t1.size(); ++index)
    {
        t1.data()[index] = std::sqrt(t2.data()[index]);
    }
}

void action5(std::vector<double>& t1, std::vector<double>& t2)
{
    #pragma omp parallel for schedule(static) num_threads(std::thread::hardware_concurrency())
    for (auto index = std::size_t{}; index < t1.size(); ++index)
    {
        t1.data()[index] = t2.data()[index] * 2.0;
    }
}

void all_actions(std::vector<double>& t1, std::vector<double>& t2)
{
    #pragma omp parallel for schedule(static) num_threads(std::thread::hardware_concurrency())
    for (auto index = std::size_t{}; index < t1.size(); ++index)
    {
        t1.data()[index] = std::sin(t2.data()[index]) * std::cos(t2.data()[index]);
        t1.data()[index] = t2.data()[index] * std::sin(t2.data()[index]);
        t1.data()[index] = t2.data()[index] * t2.data()[index];
        t1.data()[index] = std::sqrt(t2.data()[index]);
        t1.data()[index] = t2.data()[index] * 2.0;
    }
}


int main()
{
    // decide the process parameters
    const auto n = std::size_t{8000000};
    const auto test_count = std::size_t{500};
    
    // garbage data...
    auto t1 = std::vector<double>(n);
    auto t2 = std::vector<double>(n);
    
    //+/////////////////
    // perform actions one after the other
    //+/////////////////
    
    const auto sp = timer::spot_timer();
    const auto dur1 = sp.duration_in_us();
    for (auto index = std::size_t{}; index < test_count; ++index)
    {
        #pragma noinline
        action1(t1, t2);
        #pragma noinline
        action2(t1, t2);
        #pragma noinline
        action3(t1, t2);
        #pragma noinline
        action4(t1, t2);
        #pragma noinline
        action5(t1, t2);
    }
    const auto dur2 = sp.duration_in_us();
    
    //+/////////////////
    // perform all actions at once
    //+/////////////////
    const auto dur3 = sp.duration_in_us();
    for (auto index = std::size_t{}; index < test_count; ++index)
    {
        #pragma noinline
        all_actions(t1, t2);
    }
    const auto dur4 = sp.duration_in_us();
    
    const auto a = dur2 - dur1;
    const auto b = dur4 - dur3;
    if (a < b)
    {
        throw std::logic_error("negative_latency_error");
    }
    const auto fork_join_latency = (a - b) / (test_count * 4);
    
    // report
    std::cout << "Ran the program with " << omp_get_max_threads() << ", the calculated fork/join latency is: " << fork_join_latency << " us" << std::endl;
    
    return 0;
}

如您所见，我们的想法是分别执行一组操作（每个操作都在一个 OpenMP 循环中）并计算其平均持续时间，然后一起执行所有这些操作（在同一个 OpenMP 循环中) 并计算其平均持续时间。那么我们就有了两个变量的线性方程组，其中一个变量是fork/join机制的潜伏期，求解得到值

问题:

我是不是忽略了什么？
目前，我正在使用“-O0”来防止聪明的裤子编译器做它有趣的事情。我应该使用哪些编译器优化，这些是否也会对延迟本身等产生影响？
在我的 6 核 Coffee Lake 处理器上，我测得的延迟约为 850 微秒。这听起来对吗？

编辑 3

) 根据@paleonix 的建议，我在开头加入了一个热身计算，
) 为了简单起见，我减少了操作的数量，并且，
) 我已经切换到 'omp_get_wtime' 以使其易于理解。

我现在运行使用标志 -O3 来编写以下代码：

void action1(std::vector<double>& t1)
{
    #pragma omp parallel for schedule(static) num_threads(std::thread::hardware_concurrency())
    for (auto index = std::size_t{}; index < t1.size(); ++index)
    {
        t1.data()[index] = std::sin(t1.data()[index]);
    }
}

void action2(std::vector<double>& t1)
{
    #pragma omp parallel for schedule(static) num_threads(std::thread::hardware_concurrency())
    for (auto index = std::size_t{}; index < t1.size(); ++index)
    {
        t1.data()[index] =  std::cos(t1.data()[index]);
    }
}

void action3(std::vector<double>& t1)
{
    #pragma omp parallel for schedule(static) num_threads(std::thread::hardware_concurrency())
    for (auto index = std::size_t{}; index < t1.size(); ++index)
    {
        t1.data()[index] = std::atan(t1.data()[index]);
    }
}

void all_actions(std::vector<double>& t1, std::vector<double>& t2, std::vector<double>& t3)
{
    #pragma omp parallel for schedule(static) num_threads(std::thread::hardware_concurrency())
    for (auto index = std::size_t{}; index < t1.size(); ++index)
    {
        #pragma optimize("", off)
        t1.data()[index] = std::sin(t1.data()[index]);
        t2.data()[index] = std::cos(t2.data()[index]);
        t3.data()[index] = std::atan(t3.data()[index]);
        #pragma optimize("", on)
    }
}


int main()
{
    // decide the process parameters
    const auto n = std::size_t{1500000}; // 12 MB (way too big for any cache)
    const auto experiment_count = std::size_t{1000};
    
    // garbage data...
    auto t1 = std::vector<double>(n);
    auto t2 = std::vector<double>(n);
    auto t3 = std::vector<double>(n);
    auto t4 = std::vector<double>(n);
    auto t5 = std::vector<double>(n);
    auto t6 = std::vector<double>(n);
    auto t7 = std::vector<double>(n);
    auto t8 = std::vector<double>(n);
    auto t9 = std::vector<double>(n);
    
    //+/////////////////
    // warum-up, initialization of threads etc.
    //+/////////////////
    for (auto index = std::size_t{}; index < experiment_count / 10; ++index)
    {
        all_actions(t1, t2, t3);
    }
    
    //+/////////////////
    // perform actions (part A)
    //+/////////////////
    
    const auto dur1 = omp_get_wtime();
    for (auto index = std::size_t{}; index < experiment_count; ++index)
    {
        action1(t4);
        action2(t5);
        action3(t6);
    }
    const auto dur2 = omp_get_wtime();
    
    //+/////////////////
    // perform all actions at once (part B)
    //+/////////////////

    const auto dur3 = omp_get_wtime();
    #pragma nofusion
    for (auto index = std::size_t{}; index < experiment_count; ++index)
    {
        all_actions(t7, t8, t9);
    }
    const auto dur4 = omp_get_wtime();
    
    const auto a = dur2 - dur1;
    const auto b = dur4 - dur3;
    const auto fork_join_latency = (a - b) / (experiment_count * 2);
    
    // report
    std::cout << "Ran the program with " << omp_get_max_threads() << ", the calculated fork/join latency is: "
        << fork_join_latency * 1E+6 << " us" << std::endl;
    
    return 0;
}

有了这个，测得的延迟现在是 115 微秒。现在令我困惑的是，当动作改变时，这个值改变。按照我的逻辑，因为我在 A 和 B 两部分都做了同样的动作，所以实际上应该没有变化。为什么会这样？

Answer 1

这是我测量 fork-join 开销的尝试：

#include <iostream>
#include <string>

#include <omp.h>

constexpr int n_warmup = 10'000;
constexpr int n_measurement = 100'000;
constexpr int n_spins = 1'000;

void spin() {
    volatile bool flag = false;
    for (int i = 0; i < n_spins; ++i) {
        if (flag) {
            break;
        }
    }
}

void bench_fork_join(int num_threads) {
    omp_set_num_threads(num_threads);

    // create threads, warmup
    for (int i = 0; i < n_warmup; ++i) {
        #pragma omp parallel
        spin();
    }

    double const start = omp_get_wtime();
    for (int i = 0; i < n_measurement; ++i) {
        #pragma omp parallel
        spin();
    }
    double const stop = omp_get_wtime();
    double const ptime = (stop - start) * 1e6 / n_measurement;

    // warmup
    for (int i = 0; i < n_warmup; ++i) {
        spin();
    }
    double const sstart = omp_get_wtime();
    for (int i = 0; i < n_measurement; ++i) {
        spin();
    }
    double const sstop = omp_get_wtime();
    double const stime = (sstop - sstart) * 1e6 / n_measurement;

    std::cout << ptime << " us\t- " << stime << " us\t= " << ptime - stime << " us\n";
}

int main(int argc, char **argv) {
    auto const params = argc - 1;
    std::cout << "parallel\t- sequential\t= overhead\n";

    for (int j = 0; j < params; ++j) {
        auto num_threads = std::stoi(argv[1 + j]);
        std::cout << "---------------- num_threads = " << num_threads << " ----------------\n";
        bench_fork_join(num_threads);
    }

    return 0;
}

您可以使用多个不同数量的线程来调用它，线程数量不应高于您机器上的核心数量才能给出合理的结果。在我有 6 个内核并使用 gcc 11.2 编译的机器上，我得到

$ g++ -fopenmp -O3 -DNDEBUG -o bench-omp-fork-join bench-omp-fork-join.cpp
$ ./bench-omp-fork-join 6 4 2 1
parallel        - sequential    = overhead
---------------- num_threads = 6 ----------------
1.51439 us      - 0.273195 us   = 1.24119 us
---------------- num_threads = 4 ----------------
1.24683 us      - 0.276122 us   = 0.970708 us
---------------- num_threads = 2 ----------------
1.10637 us      - 0.270865 us   = 0.835501 us
---------------- num_threads = 1 ----------------
0.708679 us     - 0.269508 us   = 0.439171 us

在每一行中，第一个数字是有线程的平均值（超过 100'000 次迭代），第二个数字是没有线程的平均值。最后一个数字是前两个之间的差异，应该是 fork-join 开销的上限。

确保每一行中间列（无线程）中的数字大致相同，因为它们应该与线程数无关。如果不是，请确保计算机上没有其他东西运行 and/or 增加测量次数 and/or 预热运行。

关于将 OpenMP 替换为 MPI，请记住 MPI 仍然是多处理而非多线程。您可能会付出很多内存开销，因为进程往往比线程大得多。

编辑：

修改基准以在 volatile 标志上使用旋转而不是休眠（感谢@Jérôme Richard）。正如 Jérôme Richard 在他的回答中提到的那样，测量的开销随 n_spins 增长。将 n_spins 设置为低于 1000 并没有显着改变我的测量值，所以这就是我测量的地方。正如上面所见，测得的开销比基准测试的早期版本测得的要低得多。

休眠的不准确性是一个问题，尤其是因为人们总是会测量休眠时间最长的线程，因此会偏向于更长的时间，即使休眠时间本身会围绕输入时间对称分布。

Answer 2

TL;DR: 由于动态频率缩放，内核无法以完全相同的速度运行，并且存在大量噪声会影响执行，从而导致昂贵的同步。您的基准测试主要衡量这种同步开销。使用独特的平行部分应该可以解决这个问题。

基准测试存在相当大的缺陷。此代码实际上并未测量 OpenMP fork/join 部分的“延迟”。它衡量许多开销的组合，包括：

负载平衡和同步：拆分循环比大合并循环执行更频繁的同步（多 5 倍）。同步是昂贵的，不是因为通信开销，而是因为本质上不同步的核心之间的实际同步。实际上，线程之间的轻微 work-imbalance 会导致其他线程等待最慢线程的完成。您可能认为这不应该发生，因为静态调度，但 上下文切换 和 动态频率缩放 导致某些线程比其他线程慢。如果线程未绑定到核心，或者某些程序在计算期间由 OS 调度，则上下文切换尤为重要。 动态频率缩放（例如英特尔涡轮增压）导致一些（线程组）在工作负载方面更快，每个内核的温度和整体封装、活动核心数、预估功耗等。核心数越多，同步开销越高。请注意，此开销取决于循环所花费的时间。更多信息，请阅读以下分析。

循环拆分的性能：将 5 个循环合并成一个唯一的循环会影响生成的汇编代码（因为需要较少的指令）并且还会影响 load/store缓存（因为内存访问模式有点不同）。更不用说它在理论上会影响矢量化，尽管 ICC 不会矢量化此特定代码。话虽如此，这似乎并不是我机器上的主要实际问题，因为我无法通过运行按顺序执行程序来重现 Clang 的问题，而我可以使用多个线程。

要解决此问题，您可以使用独特的平行部分。 omp for 循环必须使用 nowait 子句，以免引入同步。或者，task-based 构造如 taskloop 和 nogroup 可以帮助实现相同的目的。在这两种情况下，您都应该注意依赖关系，因为多个 for-loop/task-loos 可以运行并行。这在您当前的代码中很好。

分析

分析由执行噪声（上下文切换、频率缩放、缓存效应、OS 中断等）引起的短同步的影响非常困难，因为它可能永远不会是最慢的线程在你的情况下同步期间（线程之间的工作是相当平衡的，但它们的速度并不完全相等）。

也就是说，如果这个假设成立，fork_join_latency 应该依赖于 n。因此，增加 n 也会增加 fork_join_latency。这是我可以在我的 6 核 i5-9600KF 处理器上使用 Clang 13 + IOMP（使用 -fopenmp -O3）：

n=   80'000    fork_join_latency<0.000001
n=  800'000    fork_join_latency=0.000036
n= 8'000'000   fork_join_latency=0.000288
n=80'000'000   fork_join_latency=0.003236

请注意，fork_join_latency 时间在实践中不是很稳定，但行为非常明显：测量的开销取决于 n。

更好的解决方案是通过测量每个线程的循环时间并累加最小和最大时间之间的差值来测量同步时间。这是一个代码示例：

double totalSyncTime = 0.0;

void action1(std::vector<double>& t1)
{
    constexpr int threadCount = 6;
    double timePerThread[threadCount] = {0};

    #pragma omp parallel
    {
        const double start = omp_get_wtime();
        #pragma omp for nowait schedule(static) //num_threads(std::thread::hardware_concurrency())
        #pragma nounroll
        for (auto index = std::size_t{}; index < t1.size(); ++index)
        {
            t1.data()[index] = std::sin(t1.data()[index]);
        }
        const double stop = omp_get_wtime();
        const double threadLoopTime = (stop - start);
        timePerThread[omp_get_thread_num()] = threadLoopTime;
    }

    const double mini = *std::min_element(timePerThread, timePerThread+threadCount);
    const double maxi = *std::max_element(timePerThread, timePerThread+threadCount);
    const double syncTime = maxi - mini;
    totalSyncTime += syncTime;
}

然后您可以按照与 fork_join_latency 相同的方式划分 totalSyncTime 并打印结果。我得到 0.000284 和 fork_join_latency=0.000398（n=8'000'000），这几乎证明了开销的主要部分是由于同步，尤其是由于线程执行速度略有不同。请注意，此开销不包括 OpenMP 并行部分末尾的隐式屏障。

Answer 3

查看我对相关问题的回答：

TLDR：我将 10k 并行循环拆分为并行区域外的 x 和内部的 10k/x。结论是，启动一个parallel region的成本基本是zip。

测量 OpenMP Fork/Join 延迟

Measuring OpenMP Fork/Join latency

c++

multithreading

openmp

performance-testing

latency

分析