OpenMP：写入与 HDF5 同步的数据

Question

我目前正在进行一个项目，其中必须使用 HDF5 创建一个大型数据集。现在，天真的实现一切都很好而且很花哨，但是很慢。慢的部分在计算中（比写入慢 10 倍），我不能再加快速度了，但也许并行化是可能的。

我想我可以使用一个简单的 #pragma omp parallel for 但出于速度原因，dataspace.write(..) 方法应该是顺序的（也许没关系）。以这张图片为例。

需要注意的是，由于维度的原因，write函数使用了一个与缓冲区大小相同的分块布局（实际上是1Mb左右）

/*
------------NAIVE IMPLEMENTATION-----------------
|T:<calc0><W0><calc1><W1><calc2><W2>............|
|-----------------------------------------------|
|----------PARALLEL IMPLEMENTATION--------------|
|-----------------------------------------------|
|T0:<calc0----><W0><calc4>.....<W4>.............|
|T1:<calc1---->....<W1><calc5->....<W5>.........|
|T2:<calc2--->.........<W2>calc6-->....<W6>.....|
|T3:<calc3----->...........<W3><calc7-->...<W7>.|
------------DIFFERENT IMPLEMENTATION-------------
i.e.: Queuesize=4
T0:.......<W0><W1><W2><W3><W4><W5><W6>..........|
T1:<calc0><calc3>.....<calc6>...................|
T2:<calc1>....<calc4>.....<calc7>...............|
T3:<calc2>........<calc5>.....<calc8>...........|


T          Thread
<calcn---> Calculation time
<Wn>       Write data n. Order *important*
.          Waiting
*/

代码示例：

#include <chrono>
#include <cmath>
#include <iostream>
#include <memory>

double calculate(float *buf, const struct options *opts) {
  // dummy function just to get a time reference
  double res = 0;
  for (size_t i = 0; i < 10000; i++)
    res += std::sin(i);
  return 1 / (1 + res);
}

struct options {
  size_t idx[6];
};

class Dataspace {
public:
  void selectHyperslab(){}; // selects region in disk space
  void write(float *buf){}; // write buf to selected disk space
};

int main() {
  size_t N = 6;
  size_t dims[6] = {4 * N, 4 * N, 4 * N, 4 * N, 4 * N, 4 * N},
         buf_offs[6] = {4, 4, 4, 4, 4, 4};
  // dims: size of each dimension, multiple of 4
  // buf_offs: size of buffer in each dimension

  // Calcuate buffer size and allocate
  // the size of the buffer is usually around 1Mb
  // and not a float but a compund datatype
  size_t buf_size = buf_offs[0];
  for (auto off : buf_offs)
    buf_size *= off;
  std::unique_ptr<float[]> buf{new float[buf_size]};

  struct options opts;        // options parameters, passed to calculation fun
  struct Dataspace dataspace; // dummy Dataspace. Supplied by HDF5

  size_t i = 0;
  size_t idx0, idx1, idx2, idx3, idx4, idx5;
  auto t_start = std::chrono::high_resolution_clock::now();
  std::cout << "[START]" << std::endl;
  for (idx0 = 0; idx0 < dims[0]; idx0 += buf_offs[0])
    for (idx1 = 0; idx1 < dims[1]; idx1 += buf_offs[1])
      for (idx2 = 0; idx2 < dims[2]; idx2 += buf_offs[2])
        for (idx3 = 0; idx3 < dims[3]; idx3 += buf_offs[3])
          for (idx4 = 0; idx4 < dims[4]; idx4 += buf_offs[4])
            for (idx5 = 0; idx5 < dims[5]; idx5 += buf_offs[5]) {
              i++;
              opts.idx[0] = idx0;
              opts.idx[1] = idx1;
              opts.idx[2] = idx2;
              opts.idx[3] = idx3;
              opts.idx[4] = idx4;
              opts.idx[5] = idx5;

              dataspace.selectHyperslab(/**/); // function from HDF5
              calculate(buf.get(), &opts);     // populate buf with data
              dataspace.write(buf.get());      // has to be sequential
            }
  std::cout << "[DONE] " << i << " calls" << std::endl;
  std::chrono::duration<double> diff =
      std::chrono::high_resolution_clock::now() - t_start;
  std::cout << "Time: " << diff.count() << std::endl;
  return 0;
}

代码应该开箱即用。

我已经快速了解了 OpenMP，但我还不能完全理解。谁能给我一个提示/工作示例？我不擅长并行化，但是带有缓冲队列的写入器线程不会工作吗？还是无论如何都使用 OpenMP 矫枉过正而 pthreads 就足够了？ 感谢您提供任何帮助，

干杯

Answer 1

您的第一个并行实现想法是迄今为止最容易实现的。制作队列和专用 I/O 线程可能性能更好，但使用 OpenMP 实现起来要困难得多。

以下是并行版本的简单示例。最重要的方面是：

1. 共享数据：确保线程间共享的任何数据都没有竞争条件。例如，每个线程必须有它自己的buf 和opts，因为它们显然是并行修改的，没有任何限制。最简单的方法是在并行区域内本地定义变量。同样循环idxn，至少对于内部循环，并且i 必须在本地定义。您不能像以前那样计算 i - 这会在每次循环迭代之间产生依赖关系并阻止并行化。
2. 将pragma omp for 工作共享应用于循环。由于每个维度的迭代次数很少，建议申请collapse。这将分配多个嵌套循环的工作。 collapse 的最佳值将为您的程序可用线程数提供足够的并行工作，但不会产生过多的开销或阻碍内部循环的单线程优化。您可能想尝试不同的值。
3. 使用critical 部分保护写入数据。一次只有一个线程将进入该部分。这很可能是正确性所必需的（取决于它在 hdf5 中的实现方式）。显然selectHyperslab 将控制write 的操作方式，因此它必须位于同一临界区。

放在一起，它可能看起来像这样：

#pragma omp parallel
{
  // define EVERYTHING that is modified locally to each thread!
  std::unique_ptr<float[]> buf{new float[buf_size]};
  struct options opts;
  // Try different values for collapse if performance is not satisfactory
  #pragma omp for collapse(3)
  for (size_t idx0 = 0; idx0 < dims[0]; idx0 += buf_offs[0])
    for (size_t idx1 = 0; idx1 < dims[1]; idx1 += buf_offs[1])
      for (size_t idx2 = 0; idx2 < dims[2]; idx2 += buf_offs[2])
        for (size_t idx3 = 0; idx3 < dims[3]; idx3 += buf_offs[3])
          for (size_t idx4 = 0; idx4 < dims[4]; idx4 += buf_offs[4])
            for (size_t idx5 = 0; idx5 < dims[5]; idx5 += buf_offs[5]) {
              size_t i = idx5 + idx4 * dims[5] + ...;
              opts.idx[0] = idx0;
              opts.idx[1] = idx1;
              opts.idx[2] = idx2;
              opts.idx[3] = idx3;
              opts.idx[4] = idx4;
              opts.idx[5] = idx5;

              calculate(buf.get(), &opts);     // populate buf with data
              #pragma omp critical
              {
                // I do assume that this function selects where/how data 
                // will be written so you *must* protected it
                // Only one thread can do this at a time.
                dataspace.selectHyperslab(/**/); // function from HDF5
                dataspace.write(buf.get());      // has to be sequential
              }
            }
}