Использование gnu FFTW 3.3.8 для вычисления одномерного FFT со сложным вводом - PullRequest
Новая
       33

Использование gnu FFTW 3.3.8 для вычисления одномерного FFT со сложным вводом

0 голосов
/ 22 июня 2019

Я использовал gnu FFTW 3.3.8 для вычисления БПФ одного измерения со сложными входами, и по некоторым причинам все работает не так, как описано ниже.

При использовании FFTW (FFTW 3.3.8 из http://www.fftw.org/) при N <= 16 код работает нормально (N - длина массива и все элементы равны единице).При N> = 22 вход и выход возвращаются в виде нулей. Мне кажется, что между двумя значениями, например, 16 

program hello
    implicit none
    integer N
    parameter (N = 20)
    double complex, dimension (N) :: in
    double complex, dimension (N) :: out
    real pi, one
    integer i
    double precision xone

    xone=1.0000000000
     pi=4.0*atan(one)


    ! generate input
    do i=1,N
         in(i)=xone
    end do

    call calc_fft(N, in, out)

    ! print output
    do i=1,N
        write(*,*)real (in(i)), real (out(i))
    end do

     ! output data into a file 
   open(1, file = 'dataM.dat', status='new')  
   do i = 1,N  
      write(1,*) in(i), out(i)   
   end do  
   close(1) 



end program hello

subroutine calc_fft(N,in,out)
    integer N
    double complex, dimension (N) :: in
    double complex, dimension (N) :: out
    integer*8 plan
    integer i

    call dfftw_plan_dft_1d(plan,N,in,out,FFTW_FORWARD,FFTW_ESTIMATE)
    call dfftw_execute_dft(plan, in, out)
    call dfftw_destroy_plan(plan)

end subroutine

gfortran testfftw.f90 -L / usr / lib64 / -lfftw3

Ввод / вывод для N = 16 

   1.0000000000000000        16.000000000000000     
   1.0000000000000000        0.0000000000000000     
   1.0000000000000000        0.0000000000000000     
   1.0000000000000000        0.0000000000000000     
   1.0000000000000000        0.0000000000000000     
   1.0000000000000000        0.0000000000000000     
   1.0000000000000000        0.0000000000000000     
   1.0000000000000000        0.0000000000000000     
   1.0000000000000000        0.0000000000000000     
   1.0000000000000000        0.0000000000000000     
   1.0000000000000000        0.0000000000000000     
   1.0000000000000000        0.0000000000000000     
   1.0000000000000000        0.0000000000000000     
   1.0000000000000000        0.0000000000000000     
   1.0000000000000000        0.0000000000000000     
   1.0000000000000000        0.0000000000000000     

 gfortran testfftw.f90 -L/usr/lib64/ -lfftw3

Input / Output for N=20
   0.0000000000000000        0.0000000000000000     
   0.0000000000000000        0.0000000000000000     
   0.0000000000000000        0.0000000000000000     
   0.0000000000000000        0.0000000000000000     
   0.0000000000000000        0.0000000000000000     
   0.0000000000000000        0.0000000000000000     
   0.0000000000000000        0.0000000000000000     
   0.0000000000000000        0.0000000000000000     
   0.0000000000000000        0.0000000000000000     
   0.0000000000000000        0.0000000000000000     
   0.0000000000000000        0.0000000000000000     
   0.0000000000000000        0.0000000000000000     
   0.0000000000000000        0.0000000000000000     
   0.0000000000000000        0.0000000000000000     
   0.0000000000000000        0.0000000000000000     
   0.0000000000000000        0.0000000000000000     
   0.0000000000000000        0.0000000000000000     
   0.0000000000000000        0.0000000000000000     
   0.0000000000000000        0.0000000000000000     
   0.0000000000000000        0.0000000000000000

Я использую Scientific Linux выпуск 7.6 (Азот) 

uname -a : 

Linux localhost.localdomain 3.10.0-957.el7.x86_64 #1 SMP Tue Oct 30 14:13:26 CDT 2018 x86_64 x86_64 x86_64 GNU/Linux

gcc-gfortran-4.8.5

Процессор Intel® Core® TM i3-3240 @ 3,40 ГГц

с 4 ГБ ОЗУ

Я также пробовал другой процессор с аналогичными эффектами

1 Ответ

1 голос
/ 23 июня 2019

Одна большая проблема кажется, что calc_fft() не имеет implicit none, поэтому применяется неявная типизация ... 

subroutine calc_fft(N,in,out)
    integer N
    double complex, dimension (N) :: in
    double complex, dimension (N) :: out
    ...

Если мы добавим implicit none как 

subroutine calc_fft(N,in,out)
    implicit none  !<--
    integer N
    double complex, dimension (N) :: in
    double complex, dimension (N) :: out
    ...

тогда gfortran выдает сообщение 

testfftw.f90:37:67:

     call dfftw_plan_dft_1d(plan,N,in,out,FFTW_FORWARD,FFTW_ESTIMATE)
                                                                   1
Error: Symbol 'fftw_estimate' at (1) has no IMPLICIT type
testfftw.f90:37:53:

     call dfftw_plan_dft_1d(plan,N,in,out,FFTW_FORWARD,FFTW_ESTIMATE)
                                                     1
Error: Symbol 'fftw_forward' at (1) has no IMPLICIT type

Здесь FFTW_FORWARD и FFTW_ESTIMATE - это некоторые параметры, которые необходимо определить с помощью заголовочного файла FFTW (в противном случае эти параметры рассматриваются как переменные real по умолчанию без implicit none!). 

subroutine calc_fft(N, in, out)
    implicit none
    include 'fftw3.f'  !<--
    integer N

Тогда мы перекомпилируем как 

gfortran testfftw.f90 -L/usr/lib64 -I/usr/include -lfftw3

и получите ожидаемый результат. (Расположение включаемых файлов может различаться в зависимости от машин / ОС, и ScientificLinux7, кажется, имеет их в /usr/include. Пожалуйста, посмотрите эту страницу и, при необходимости, установите их как # yum install fftw-devel. Также чтобы получить лучшую производительность, различные пакеты FFTW могут быть лучше, чем полученные из yum, но это другая история ...)

...