首页 > 解决方案 > 使用 gnu FFTW 3.3.8 计算具有复杂输入的一维 FFT

问题描述

我使用 gnu FFTW 3.3.8 来计算具有复杂输入的一维 FFT,由于某些原因,如下所述,事情无法正常工作。

当使用 N<=16 的 FFTW(来自http://www.fftw.org/的 FFTW 3.3.8 )时,代码运行良好(N 是数组长度,所有元素都等于 1)。对于 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(R) Core(TM) i3-3240 CPU @ 3.40GHz

带 4GB 内存

我还尝试了另一个具有类似效果的 CPU

标签: fortrangfortranfftw

解决方案


一个大问题似乎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_FORWARDFFTW_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-develyum


推荐阅读