milesdavis
milesdavis

Reputation: 1

Using gnu FFTW 3.3.8 to compute one dimensional FFT with complex input

I used gnu FFTW 3.3.8 to calculate one dimension FFT with complex inputs and for some reasons things are not working as described below.

When using FFTW ( FFTW 3.3.8 from http://www.fftw.org/ ) with N<=16 the code runs fine (N is the array length and all elements are equal to one). For N>=22 both the input and output are returned as zeros.It seems to me that in between the two values e.g. 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

Input / Output for 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     

I am using Scientific Linux release 7.6 (Nitrogen)

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

with 4gb of ram

I also tried on another CPU with similar effects

Upvotes: 0

Views: 221

Answers (1)

roygvib
roygvib

Reputation: 7395

One big issue seems that calc_fft() does not have implicit none, so implicit typing applies...

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

If we add implicit none as

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

then gfortran gives the message

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

Here, FFTW_FORWARD and FFTW_ESTIMATE are some parameters that need to be defined via a header file of FFTW (otherwise, those parameters are regarded as default real variables without implicit none!).

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

Then we recompile as

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

and get the expected result. (The location of the include files may vary depending on machines/OS, and ScientificLinux7 seems to have them in /usr/include. Please see this page and, if necessary, install them as # yum install fftw-devel. Also, to get the best performance, different FFTW packages might be better than that obtained from yum, but it is a different story...)

Upvotes: 1

Related Questions