在具有多线程的 OOP 中使用 FFTW

问题描述

考虑一个二维复数值数组的容器

#include <vector>
#include <array>

struct Array2D {
  typedef std::array<double,2> complex;

  int X,Y;
  std::vector<complex> values;
};

以及使用 FFTW 计算傅立叶变换的函数

#include <fftw3.h>

void FourierTransform(const Array2D& source,Array2D *dest) {
  fftw_plan p = fftw_plan_dft_2d(source.X,source.Y,(fftw_complex*)source.values.data(),(fftw_complex*)dest->values.data(),FFTW_FORWARD,FFTW_ESTIMATE);
  fftw_execute(p);
  fftw_destroy_plan(p);
}

这不是最佳实现,因为它在每次计算傅立叶变换时都会创建和销毁计划,但它方便且易于使用。但是,创建和销毁 FFTW 计划不是线程安全的(请参阅 here),因此不应从不同线程同时调用函数

问题:

  • 我的应用程序中有许多不同的数组,我需要计算傅立叶变换和
  • 在主循环与 #pragma omp parallel for 的并行化和代码中需要计算傅立叶变换的地方之间有几个函数调用

#pragma omp parallel for 之前初始化所有这些数组以及 fftw 计划,然后将它们作为参数传递给后续函数会使代码很多不那么清晰。有没有另一种方法可以让我以线程安全的方式封装和隐藏 FFTW 函数调用

解决方法

FFTW 文档的链接提出了一个答案,即确保对 FFTW 的不安全调用被序列化。它建议使用信号量来做到这一点,但由于您使用的是 OpenMP,您可以使用 #pragma omp critical 来做到这一点,就像这样

    void FourierTransform(const Array2D& source,Array2D *dest) {
      fftw_plan p;
      #pragma omp critical (FFTW)
      {
        p = fftw_plan_dft_2d(source.X,source.Y,(fftw_complex*)source.values.data(),(fftw_complex*)dest->values.data(),FFTW_FORWARD,FFTW_ESTIMATE);
      }
      fftw_execute(p);
      #pragma omp critical (FFTW)
      {
         fftw_destroy_plan(p);
      }
    }

我在那里使用了一个命名的临界区(名称为 FFTW),以便这些临界区不会降低与代码中任何其他临界区的并行度,而是相互保护。

更一般地说,提供相关函数的包装的、线程安全的版本并调用它们可能会很好,而不是未包装的。那将是类似

static fftw_plan_p TS_fftw_plan_dft_2d(... arguments ... ) {
    fftw_plan p;
    #pragma omp critical (FFTW)
    {
        p = fftw_plan_dft_2d(... arguments ...);
    }
    return p;
}

然后您只需在其余代码中调用 TS_fftw_foo 而不是 fftw_foo

如果你想变得极端,你可以将所有这些包装器放入一个标题中(也将它们标记为内联),然后还有宏来将相关的 FFTW 函数扩展到 TS_fftw 函数中。然后包含您的标题就可以完成这项工作。 不过,这更加模糊和极端。只需包装您需要的函数可能会更清楚地了解正在发生的事情。

显然所有这些代码都未经测试和编译,但我希望你能明白。