nano_exit

基礎的なことこそ、簡単な例が必要だと思うのです。

FFTWをfortranで使う際の自分用のメモ

FFTWをfortranで使う。

本家の説明:
FFTW 3.3.8: Fortran Examples

参考にしたサイト:
fftwの使い方:腰も砕けよ 膝も折れよ:So-net blog
FortranのFFTWの使い方 - Qiita
Hiroyuki R. Takahashi
Fumio KUSUHARA -Fortran90で構造解析-

program fftw3_test
implicit none

integer,       parameter :: WP = selected_real_kind( 8 )
real( WP ),    parameter :: PI = 3.141592653589793

integer :: i1

! setting parameters and variables
integer,     parameter :: N = 2**10
real( WP ) , parameter :: L = 100.
real( WP )    :: dx
real( WP )    :: x( N ), k( N )
complex( WP ) :: diff

! data
complex( WP ) :: input0( N ), input( N ), output( N )

! plan
integer( WP ) :: plan

include 'fftw3.f'

dx = L / N

! test: Gaussian function
do i1 = 1, N
  x( i1 )     = - L/2 + ( i1 - 1 ) * dx
  input( i1 ) = exp( - x( i1 )**2 )
enddo

input0( : ) = input( : )

call dfftw_plan_dft_1d( plan, N, input0, output,    &
  &                     FFTW_FORWARD, FFTW_ESTIMATE )
call dfftw_execute_dft( plan, input0, output )
call dfftw_destroy_plan( plan )

do i1 = 1, N/2
  k( i1             ) =   i1 - 1
  k( N - ( i1 - 1 ) ) = - i1
enddo
k( : ) = ( 2 * PI / L ) *  k( : ) 

open( 10, file = 'fftw3_forward.out', status = 'unknown' )
do i1 = 1, N
  write( 10, '(4e15.5)' ) k( i1 ), &
    &  real( output( i1 ) ), aimag( output( i1 ) ), abs( output( i1 ) )
enddo
close( 10 )

call dfftw_plan_dft_1d( plan, N, output, input,     &
  &                     FFTW_BACKWARD, FFTW_ESTIMATE )
call dfftw_execute_dft( plan, output, input )
call dfftw_destroy_plan( plan )
input( : ) = input( : ) / N

open( 10, file = 'fftw3_backward.out', status = 'unknown' )
do i1 = 1, N
  diff = input( i1 ) - input0( i1 )
  write( 10, '(4e15.5)' ) x( i1 ), &
    &  real( diff ), aimag( diff ), abs( diff )
enddo
close( 10 )

end program fftw3_test

フーリエ変換後の配列は、0, 1, 2, ..., N/2, -(N/2 - 1), -(N/2 - 2), ... -1 となっていることに注意。
この配列のまま逆フーリエ変換にかけると元に戻る。

compile:

gfortran -o fftw3_test.out fftw3_test.f90 -I/usr/local/include -L/usr/local/lib -lfftw3