openmp - FFTW3 - 就地并行化一维复数 fft 很慢
问题描述
所以我正在研究并行化一维 FFT。作为第一项任务,我在Intel(R) Xeon(R) CPU E5-2620 v3 @ 2.40GHz(有 16 个内核)上执行了 FFTW3 库的基准测试。我刚刚做了一个基本的一维复数 FFT,用 OpenMP 作为我的线程库。我使用以下命令在 ICC 上编译:
icc -Wall -Werror
-I/.../mkl/include -I/apps/intel/linux/mkl/include/fftw
fftw3_dft.c
-L/.../intel/linux/mkl/.../intel64 -lmkl_rt
-L/.../intel/.../linux/mkl/../compiler/lib/intel64
-L/apps/intel/.../clinux/mkl/../tbb/lib/intel64/gcc4.4
-liomp5 -lm -lpthread -ldl
-o fftw3_dft.out
我计算了不同问题大小的加速指标。我无法解释这个情节
- 对于 2^21 和 2^24 之间的问题大小,为什么使用 2 个处理器没有加速?(即使 4,8 和 16 线程有一些加速)
- 当问题大小变得大于 2^27 时,为什么会突然增加加速?
代码
#include <stdio.h>
#include <math.h>
#include <stdlib.h>
#include <float.h>
#include "fftw3.h"
#include "mkl.h"
/* Compute (K*L)%M accurately */
static double moda(int K, int L, int M)
{
return (double)(((long long)K * L) % M);
}
/* Initialize array x[N] with harmonic H */
static void init(fftw_complex *x, int N, int H)
{
double TWOPI = 6.2831853071795864769, phase;
int n;
for (n = 0; n < N; n++)
{
phase = moda(n,H,N) / N;
x[n][0] = cos( TWOPI * phase ) / N;
x[n][1] = sin( TWOPI * phase ) / N;
}
}
int main(int argc, char *argv[]) {
if(argc < 3) {
printf("Error : give args\n");
return 0;
}
int N = atoi(argv[1]);
int p = atoi(argv[2]);
int H = -N/2;
fftw_plan forward_plan = 0, backward_plan = 0;
fftw_complex *x = 0;
int status = 0;
fftw_init_threads();
fftw_plan_with_nthreads(p);
x = fftw_malloc(sizeof(fftw_complex)*N);
forward_plan = fftw_plan_dft(1, &N, x, x, FFTW_FORWARD, FFTW_ESTIMATE);
init(x, N, H);
double start_time = dsecnd();
/*--------------ALG STARTS HERE --------------------------*/
fftw_execute(forward_plan);
/*--------------ALG ENDS HERE --------------------------*/
double end_time = dsecnd();
printf(LI", %d, %lf\n", N, p, end_time - start_time);
fftw_cleanup_threads()
fftw_destroy_plan(forward_plan);
fftw_free(x);
}
解决方案
即使在 n=2^14 的情况下,我也确实获得了 2 个内核的加速,之后它始终保持在 1.5 以上。请记住多次运行代码并丢弃用于启动的第一部分。现代内核需要一些时间才能全速运行。
#include <stdio.h>
#include <math.h>
#include <stdlib.h>
#include <float.h>
#include "fftw3.h"
#include "omp.h"
/* Compute (K*L)%M accurately */
static double moda(int K, int L, int M)
{
return (double)(((long long)K * L) % M);
}
/* Initialize array x[N] with harmonic H */
static void init(fftw_complex *x, int N, int H)
{
double TWOPI = 6.2831853071795864769, phase;
int n;
for (n = 0; n < N; n++)
{
phase = moda(n,H,N) / N;
x[n][0] = cos( TWOPI * phase ) / N;
x[n][1] = sin( TWOPI * phase ) / N;
}
}
int main(int argc, char *argv[]) {
if(argc < 2) {
printf("Error : give args\n");
return 0;
}
int max_pow = atoi(argv[1]);
int p = 1;
#pragma omp parallel
{
#pragma omp single
{
p = omp_get_num_threads();
}
}
printf("%i\n", p);
fftw_plan forward_plan = 0;
fftw_complex *x = 0;
fftw_init_threads();
fftw_plan_with_nthreads(p);
for(int iter=1;iter<=2;iter++){
//throw away the first round, a couple of seconds is enough
for(int pw=12;pw<=max_pow;pw++){
int N = pow(2, pw);
int H = -N/2;
x = fftw_malloc(sizeof(fftw_complex)*N);
forward_plan = fftw_plan_dft(1, &N, x, x, FFTW_FORWARD, FFTW_MEASURE);
init(x, N, H);
double start_time = omp_get_wtime();
/*--------------ALG STARTS HERE --------------------------*/
for(int i=1;i<=5;i++){fftw_execute(forward_plan);}
/*--------------ALG ENDS HERE --------------------------*/
double end_time = omp_get_wtime();
printf("%i %lf\n", pw, (end_time - start_time)/5);
fftw_destroy_plan(forward_plan);
fftw_free(x);
}
}
return 0;
}
和
> gfortran -fopenmp fftw1d.c -lfftw3 -lfftw3_omp
> OMP_NUM_THREADS=2 ./a.out 24
结果与-O3
已经编译的库相同。
在四核 Intel(R) Core(TM) i7-3770 CPU @ 3.40GHz 上测试
推荐阅读
- c# - 是否可以在 IEntityTypeConfiguration 中引用相关实体
配置方法? - php - Phpspreadsheet 转换后的文件应具有包含相同值的附加新列
- c++ - 如何找到当前硬件线程的 L3 缓存索引和 NUMA 节点索引
- python - 使用字典值创建新列
- gstreamer - 使用 gstreamer 获得确定性时间
- c# - 列表框项目初始化,但一旦选择了一个,我不知道如何重置没有选择任何项目的列表框
- javascript - 如何运行每 30 秒调用一次 axios 的函数
- java - 使用毕加索从 Firebase 向 CardView 显示用户图像不起作用
- spartacus-storefront - 本地 Commerce-Spartacus 构建失败 - 无法加载 shrinkwap
- arangodb - How add new key with update in Arangodb