使用OpenMP创建FFTW计划

Question

我正在尝试并行执行几个FFT。我正在使用FFTW和OpenMP。每个FFT都是不同的，所以我不依赖FFTW的内置多线程（我知道它使用OpenMP）。

int m;
// assume:
// int numberOfColumns = 100;
// int numberOfRows = 100;
#pragma omp parallel for default(none) private(m) shared(numberOfColumns, numberOfRows)//  num_threads(4)
    for(m = 0; m < 36; m++){
        // create pointers
        double          *inputTest;
        fftw_complex    *outputTest;
        fftw_plan       testPlan;
        // preallocate vectors for FFTW
         outputTest = (fftw_complex*)fftw_malloc(sizeof(fftw_complex)*numberOfRows*numberOfColumns);
         inputTest  = (double *)fftw_malloc(sizeof(double)*numberOfRows*numberOfColumns);
         // confirm that preallocation worked
         if (inputTest == NULL || outputTest == NULL){
             logger_.log_error("tt FFTW memory not allocated on m = %i", m);
         }
         // EDIT: insert data into inputTest
         inputTest = someDataSpecificToThisIteration(m); // same size for all m
        // create FFTW plan
        #pragma omp critical (make_plan)
        {
            testPlan = fftw_plan_dft_r2c_2d(numberOfRows, numberOfColumns, inputTest, outputTest, FFTW_ESTIMATE);
        }
         // confirm that plan was created correctly
         if (testPlan == NULL){
             logger_.log_error("tt failed to create plan on m = %i", m);
         }
        // execute plan
         fftw_execute(testPlan);
        // clean up
         fftw_free(inputTest);
         fftw_free(outputTest);
         fftw_destroy_plan(testPlan);
    }// end parallelized for loop

这一切都很好。但是，如果我在创建计划时删除关键结构（fftw_plan_dft_r2c_2d），我的代码就会失败。有人能解释一下原因吗？fftw_plan_dft_r2c_2d并不是真正的"孤儿"，对吧？是因为两个线程可能都试图同时命中numberOfRows或numberOfColumns

Answer 1

这几乎都写在关于线程安全的FFTW文档中：

但是必须小心，因为计划程序在调用和计划之间共享数据（例如智慧和三角表）。

结果是，FFTW中唯一的线程安全（可重入）例程是fftw_execute（及其新的数组变体）。所有其他例程（例如计划器）一次只能从一个线程调用。因此，例如，您可以将信号量锁封装在对计划器的任何调用周围；更简单的是，您可以从一个线程创建所有计划。我们认为这不应该是一个重要的限制（FFTW是为唯一对性能敏感的代码是转换的实际执行的情况而设计的），并且计划之间共享数据的好处是巨大的。

在FFT的典型应用中，很少构建计划，因此是否必须同步创建它们并不重要。在您的情况下，您不需要在每次迭代中创建新的计划，除非数据的维度发生变化。您宁愿执行以下操作：

#pragma omp parallel default(none) private(m) shared(numberOfColumns, numberOfRows)
{
   // create pointers
   double          *inputTest;
   fftw_complex    *outputTest;
   fftw_plan       testPlan;
   // preallocate vectors for FFTW
   outputTest = (fftw_complex*)fftw_malloc(sizeof(fftw_complex)*numberOfRows*numberOfColumns);
   inputTest  = (double *)fftw_malloc(sizeof(double)*numberOfRows*numberOfColumns);
   // confirm that preallocation worked
   if (inputTest == NULL || outputTest == NULL){
      logger_.log_error("tt FFTW memory not allocated on m = %i", m);
   }
   // create FFTW plan
   #pragma omp critical (make_plan)
   testPlan = fftw_plan_dft_r2c_2d(numberOfRows, numberOfColumns, inputTest, outputTest, FFTW_ESTIMATE);
   #pragma omp for
   for (m = 0; m < 36; m++) {
      // execute plan
      fftw_execute(testPlan);
   }
   // clean up
   fftw_free(inputTest);
   fftw_free(outputTest);
   fftw_destroy_plan(testPlan);
}

现在，计划在每个线程中只创建一次，串行化开销将随着fftw_execute()的每次执行而减少。如果在NUMA系统（例如，多插槽AMD64或Intel（post-）Nehalem系统）上运行，则应启用线程绑定以实现最大性能。