Главная » C++ » Как использовать интерфейс fftw Guru

Как использовать интерфейс fftw Guru

Я использовал, чтобы использовать fftw_plan_dft для многомерного преобразования Фурье. 

fftw_plan fftw_plan_dft(int rank, const int *n, fftw_complex *in,
fftw_complex *out, int sign, unsigned flags);

Теперь я хочу передать 64-битное целое число в fftw, похоже, мне нужно использовать интерфейс гуру fftw.

 fftw_plan fftw_plan_guru64_dft(
int rank, const fftw_iodim64 *dims,
int howmany_rank, const fftw_iodim64 *howmany_dims,
fftw_complex *in, fftw_complex *out,
int sign, unsigned flags);

Но я не понимаю, что это howmany_rank а также howmany_dims имею в виду. Руководство по fftw_plan_guru_dft говорит:

Эти две функции планируют многомерный ДПФ со сложными данными для чередующегося и разделенного форматов соответственно. Размеры преобразования задаются как (rank, dims) над многомерным вектором (loop) измерений (howmany_rank, howmany_dims). dims и howmany_dims должны указывать на массивы fftw_iodim ранга длины и howmany_rank соответственно.

Я знаю, знаю, что означает «многомерный вектор (цикл) измерений (howmany_rank, howmany_dims)». Можете ли вы дать мне пример или объяснить, как использовать этот интерфейс гуру?

3

Решение

Если размеры и размеры ваших многомерных массивов больше 2 ^ 32, 64-битный интерфейс гуру становится полезным.

Прототип функции создания сложных в сложные DTF:

fftw_plan fftw_plan_guru64_dft(
int rank, const fftw_iodim64 *dims,
int howmany_rank, const fftw_iodim64 *howmany_dims,
fftw_complex *in, fftw_complex *out,
int sign, unsigned flags);

где:

  • rank является рангом преобразования FFTW, которое должно быть выполнено, то есть количеством измерений.
  • dims массив размеров rank, Для каждого измерения i, dims[i].n размер линии, dims[i].is шаг между строками входного массива и dims[i].os шаг между строками выходного массива Например, если массив является непрерывным в памяти, то документация интерфейса гуру предлагает использовать рецидив dims[i].is = n[i+1] * dims[i+1].is,
    Количество выполняемых преобразований и смещения между начальными точками определяются как howmany_rank а также howmany_dims,
  • howmany_rank определяет количество преобразований с определенными смещениями.
  • howmany_dims массив размеров howmany_rank, Для каждого преобразования i, howmany_dims[i].n число преобразований, которые должны быть вычислены, каждое из которых имеет смещение между входами howmany_dims[i].is и смещение между выходом howmany_dims[i].os,

Следующий код вызывает fftw_plan_guru64_dft() так что он выполняет то же самое, что и fftw_plan_dft_3d(), Это может быть скомпилировано gcc main.c -o main -lfftw3 -lm -Wall:

#include<stdlib.h>
#include<complex.h>
#include<math.h>
#include<fftw3.h>

int main(void){

fftw_plan p;
unsigned long int N = 10;
unsigned long int M = 12;
unsigned long int P = 14;
fftw_complex *in=fftw_malloc(N*M*P*sizeof(fftw_complex));
if(in==NULL){fprintf(stderr,"malloc failed\n");exit(1);}
fftw_complex *out=fftw_malloc(N*M*P*sizeof(fftw_complex));
if(out==NULL){fprintf(stderr,"malloc failed\n");exit(1);}
unsigned int i,j,k;

int rank=3;
fftw_iodim64 *dims=malloc(rank*sizeof(fftw_iodim64));
if(dims==NULL){fprintf(stderr,"malloc failed\n");exit(1);}
dims[0].n=N;
dims[0].is=P*M;
dims[0].os=P*M;
dims[1].n=M;
dims[1].is=P;
dims[1].os=P;
dims[2].n=P;
dims[2].is=1;
dims[2].os=1;

int howmany_rank=1;
fftw_iodim64 *howmany_dims=malloc(howmany_rank*sizeof(fftw_iodim64));
if(howmany_dims==NULL){fprintf(stderr,"malloc failed\n");exit(1);}
howmany_dims[0].n=1;
howmany_dims[0].is=1;
howmany_dims[0].os=1;

printf("sizeof fftw complex %ld\n",sizeof(fftw_complex));
printf("sizeof fftw_iodim64 %ld\n",sizeof(fftw_iodim64));
printf("creating the plan\n");
p=fftw_plan_guru64_dft(rank, dims,howmany_rank, howmany_dims,in, out,FFTW_FORWARD, FFTW_ESTIMATE);
if (p==NULL){fprintf(stderr,"plan creation failed\n");exit(1);}
printf("created the plan\n");

for(i=0;i<N;i++){
for(j=0;j<M;j++){
for(k=0;k<P;k++){
//printf("ijk\n");
in[(i*M+j)*P+k]=30.+12.*sin(2*3.1415926535*i/((double)N))*sin(2*3.1415926535*j/((double)M))*sin(2*3.1415926535*k/((double)P))*I;
}
}
}

fftw_execute(p);

for (i = 0; i < N; i++){
for (j = 0; j < M; j++){
for (k = 0; k < P; k++){
printf("result: %d %d %d %g %gI\n", i,j,k, creal(out[(i*M+j)*P+k]), cimag(out[(i*M+j)*P+k]));
}
}
}fftw_destroy_plan(p);
fftw_free(in);
fftw_free(out);

free(dims);
free(howmany_dims);

return(0);
}

Например, интерфейс гуру может использоваться для вычисления ДПФ сложного трехмерного электрического поля. В каждой точке сетки электрическое поле представляет собой вектор размером 3. Следовательно, я могу хранить электрическое поле в виде 4-мерной матрицы, причем последнее измерение определяет компонент вектора. Наконец, интерфейс гуру можно использовать для одновременного выполнения трех трехмерных ДПФ:

#include<stdlib.h>
#include<complex.h>
#include<math.h>
#include<fftw3.h>

int main(void){

fftw_plan p;
unsigned long int N = 10;
unsigned long int M = 12;
unsigned long int P = 14;
unsigned long int DOF = 3;
fftw_complex *in=fftw_malloc(N*M*P*DOF*sizeof(fftw_complex));
if(in==NULL){fprintf(stderr,"malloc failed\n");exit(1);}
fftw_complex *out=fftw_malloc(N*M*P*DOF*sizeof(fftw_complex));
if(out==NULL){fprintf(stderr,"malloc failed\n");exit(1);}
unsigned int i,j,k;

int rank=3;
fftw_iodim64 *dims=malloc(rank*sizeof(fftw_iodim64));
if(dims==NULL){fprintf(stderr,"malloc failed\n");exit(1);}
dims[0].n=N;
dims[0].is=P*M*DOF;
dims[0].os=P*M*DOF;
dims[1].n=M;
dims[1].is=P*DOF;
dims[1].os=P*DOF;
dims[2].n=P;
dims[2].is=DOF;
dims[2].os=DOF;

int howmany_rank=1;
fftw_iodim64 *howmany_dims=malloc(howmany_rank*sizeof(fftw_iodim64));
if(howmany_dims==NULL){fprintf(stderr,"malloc failed\n");exit(1);}
howmany_dims[0].n=3;
howmany_dims[0].is=1;
howmany_dims[0].os=1;

printf("sizeof fftw complex %ld\n",sizeof(fftw_complex));
printf("sizeof fftw_iodim64 %ld\n",sizeof(fftw_iodim64));
printf("creating the plan\n");
p=fftw_plan_guru64_dft(rank, dims,howmany_rank, howmany_dims,in, out,FFTW_FORWARD, FFTW_ESTIMATE);
if (p==NULL){fprintf(stderr,"plan creation failed\n");exit(1);}
printf("created the plan\n");

for(i=0;i<N;i++){
for(j=0;j<M;j++){
for(k=0;k<P;k++){
//printf("ijk\n");
in[((i*M+j)*P+k)*3]=30.+12.*sin(2*3.1415926535*i/((double)N))*sin(2*3.1415926535*j/((double)M))*sin(2*3.1415926535*k/((double)P))*I;
in[((i*M+j)*P+k)*3+1]=42.0;
in[((i*M+j)*P+k)*3+2]=1.0;
}
}
}

fftw_execute(p);

for (i = 0; i < N; i++){
for (j = 0; j < M; j++){
for (k = 0; k < P; k++){
printf("result: %d %d %d || %g %gI | %g %gI | %g %gI\n", i,j,k, creal(out[((i*M+j)*P+k)*3]), cimag(out[((i*M+j)*P+k)*3]),creal(out[((i*M+j)*P+k)*3+1]), cimag(out[((i*M+j)*P+k)*3+1]),creal(out[((i*M+j)*P+k)*3+2]), cimag(out[((i*M+j)*P+k)*3+2]));
}
}
}fftw_destroy_plan(p);
fftw_free(in);
fftw_free(out);

free(dims);
free(howmany_dims);

return(0);
}
1

Другие решения

Других решений пока нет …

Источник
По вопросам рекламы [email protected]
Web-Answers © 2024 Наверх