t-nissie/00FFTW.md

## 00FFTW.md

      
    Raw
  

              00FFTW.md
            
          
    FFTW3で自分のコードレット (codelet) を作ってみる

この文章はFFTW3の使い方を説明するものではない。
必要なもの


gcc
gfortran
OCaml（yum install ocamlとかでインストールできる。）
autoconf
automake（apt-get install autoconf automake libtoolとかでインストールできる。）
libtool
fig2dev（MacPortsならport install transfigでインストールできる。）
w3m（Firefoxなどのブラウザでも可。wgetやcurlではリダイレクトのせいでうまくいかないようだった。）

準備

$ git clone https://github.com/FFTW/fftw3.git fftw3
$ cd fftw3
$ git branch -a
* master
  remotes/origin/HEAD -> origin/master
  remotes/origin/experimental-simd
  remotes/origin/master
$ git checkout -b mybranch origin/master # そのorigin/masterからmybranchをブランチ
$ git branch -a
  master
* mybranch
  remotes/origin/HEAD -> origin/master
  remotes/origin/experimental-simd
  remotes/origin/master
$ ./bootstrap.sh --enable-openmp --enable-sse2 --enable-avx2
$ make -j4

CPUとアセンブラが対応していれば、version 3.3.5から有効になるであろうオプション--enable-avx2が使える。
--enable-fmaは--enable-avx2の場合は自動的に有効になる？
使ってみる

倍精度の実数⇆複素数の３次元FFTで使ってみる。上の『準備』に続けて、
$ git clone https://gist.github.com/bf784a5bc070a6f60a76.git SoA
$ echo SoA >> .gitignore
$ echo '*dft/simd/avx*-128-fma/*.c' >> .gitignore
$ git status
$ cd SoA
$ w3m https://sourceforge.net/p/loto/code/HEAD/tree/feram/branches/z_star_r/src/feram_fftw_SoA_out.F?format=raw
$ w3m https://sourceforge.net/p/loto/code/HEAD/tree/feram/branches/z_star_r/src/feram_fftw_wisdom_timing.F?format=raw
$ make -j3
$ for len in 1 2 3 4 5 6; do OMP_NUM_THREADS=2 ./feram_fftw_SoA_out_len 100 64 64 64 1 $len FFTW_PATIENT 2> /dev/null; done

はじめにこのgistをcloneし、次にw3mコマンドでferamのリポジトリ
https://sourceforge.net/p/loto/code/HEAD/tree/feram/branches/z_star_r/src/ から
feram_fftw_SoA_out.Fとferam_fftw_wisdom_timing.Fとをダウンロードしている。
最後の実行結果は2.5 GHz (i5-3210M) dual-core Intel Core i5 Ivy Bridge processorでは
     100  1   64   64   64    1     262144 out   2   0.175  13.5
     100  2   64   64   64    1     262144 out   2   0.534   8.8
     100  3   64   64   64    1     262144 out   2   0.729   9.7
     100  4   64   64   64    1     262144 out   2   0.971   9.7
     100  5   64   64   64    1     262144 out   2   1.205   9.8
     100  6   64   64   64    1     262144 out   2   1.468   9.6

となる。最後がGFLOPS値。
このようにSoA (Structure of Array) でAarryの数を１より大きくすると遅くなるので、
FFTW3の内部の変更で全部が最高の値 (13.5) になるようにしたい。
簡単なはずだが、敷居が高い。
perfでプロファイリング

americium02でperfでプロファイリングした。
よいwisdomを作る

何回かferam_fftw_SoA_out_lenを実行して、最速になるwisdom_SoAをwisdomとして残す。
$ mkdir perf_len1_8k_96x96x96_1
$ cd    perf_len1_8k_96x96x96_1
$ numactl --cpubind 1 ../feram_fftw_SoA_out_len 1000 96 96 96 1 1 FFTW_PATIENT
$ numactl --cpubind 1 ../feram_fftw_SoA_out_len 1000 96 96 96 1 1 FFTW_PATIENT
$ mv wisdom_SoA wisdom
$ numactl --cpubind 1 ../feram_fftw_SoA_out_len 1000 96 96 96 1 1 FFTW_PATIENT
$ sh ../../tools/fftw-wisdom-to-conf < wisdom   # I do not know if it helps you or not.

perfでプロファイリング

10秒前後になるようにイテレーションの回数を調節してperfでプロファイリングする。
$ numactl --cpubind 1 perf record ../feram_fftw_SoA_out_len 8000 96 96 96 1 1 FFTW_PATIENT
feram_fftw_SoA_out_len.F:30: Successfully imported FFTW wisdom in current directory.
feram_fftw_SoA_out_len.F:63: flags = FFTW_PATIENT
    8000  1   96   96   96    1     884736 out  12   8.312  84.1
[ perf record: Woken up 59 times to write data ]
[ perf record: Captured and wrote 15.258 MB perf.data (~666622 samples) ]
$ perf report > perf.report

なお、kernel云々のperfがエラーを出したら、
# echo 0 > /proc/sys/kernel/kptr_restrict

とする。

  
## feram_fftw_SoA_out_len.F
! feram_fftw_SoA_out_len.F -*-f90-*-
! Time-stamp: <2016-01-04 14:43:04 t-nissie>
! Author: Takeshi NISHIMATSU
!!
#if defined HAVE_CONFIG_H
#  include "config.h"
#endif
program feram_fftw_SoA_out_len
  use, intrinsic :: iso_c_binding
  implicit none
  include 'fftw3.f03'
  integer :: N_TIMES = 100
  integer :: Lx = 64
  integer :: Ly = 64
  integer :: Lz = 64
  integer :: padding_y = 1
  integer :: length_structure = 1
  integer :: flags = FFTW_MEASURE
  character(len=100) :: str = 'FFTW_MEASURE'
  integer :: OMP_GET_MAX_THREADS
  integer :: narg, ireturn
  type(C_PTR) :: plan_r2c_1_out = c_null_ptr
  type(C_PTR) :: plan_c2r_1_out = c_null_ptr

  ireturn = fftw_init_threads()
  call fftw_plan_with_nthreads(OMP_GET_MAX_THREADS())

  ireturn = fftw_import_wisdom_from_filename(C_CHAR_"wisdom"//C_NULL_CHAR)
  if (ireturn.eq.1) then
     write(0,'(a,a,i2,a)') __FILE__, ':', __LINE__, ': Successfully imported FFTW wisdom in current directory.'
  else
     write(0,'(a,a,i2,a)') __FILE__, ':', __LINE__, ': Failed to import FFTW wisdom in current directory.'
  end if

  narg = command_argument_count()
  if (narg.eq.0) then
     ! default values
  else if (narg.eq.7) then
     call get_command_argument(1,str); read(str,*) N_TIMES
     call get_command_argument(2,str); read(str,*) Lx
     call get_command_argument(3,str); read(str,*) Ly
     call get_command_argument(4,str); read(str,*) Lz
     call get_command_argument(5,str); read(str,*) padding_y
     call get_command_argument(6,str); read(str,*) length_structure
     call get_command_argument(7,str)
     select case(str)
     case(     'FFTW_ESTIMATE')
        flags = FFTW_ESTIMATE
     case(     'FFTW_PATIENT')
        flags = FFTW_PATIENT
     case(     'FFTW_MEASURE')
        flags = FFTW_MEASURE
     case(     'FFTW_EXHAUSTIVE')
        flags = FFTW_EXHAUSTIVE
     case default
        write(0,   '(a,a,i2,a)') __FILE__, ':', __LINE__, ': No such FFTW flag.'
        stop 2
     end select
  else
     write(0,   '(a,a,i2,a)') __FILE__, ':', __LINE__, ': Illegal number of arguments.'
     stop 1
  end if
  write(0,'(a,a,i2,a)') __FILE__, ':', __LINE__, ': flags = '//trim(str)

  call feram_fftw_SoA_out(length_structure, N_TIMES, Lx, Ly, Lz, padding_y, flags, &
       &                                                        plan_r2c_1_out, &
       &                                                        plan_c2r_1_out)

#ifndef WITH_FFT_MKL
  ireturn = fftw_export_wisdom_to_filename(C_CHAR_"wisdom_SoA"//C_NULL_CHAR)
#endif
  call fftw_destroy_plan(plan_r2c_1_out)
  call fftw_destroy_plan(plan_c2r_1_out)
  call fftw_cleanup_threads()
end program feram_fftw_SoA_out_len

## Makefile
# -*-Makefile-*- for feram_fftw_SoA_out_len
##
FC=gfortran
FFLAGS=-Wall -Ofast -fopenmp -ffree-form
CPPFLAGS=-I../api
LDFLAGS=../.libs/libfftw3.a ../threads/.libs/libfftw3_omp.a
feram_fftw_SoA_out_len: feram_fftw_SoA_out_len.o feram_fftw_SoA_out.o feram_fftw_wisdom_timing.o
	$(FC) $(FFLAGS) -o $@ $^ $(LDFLAGS)
clean:
	rm -f feram_fftw_*.o feram_fftw_SoA_out_len
	! feram_fftw_SoA_out_len.F --f90--
	! Time-stamp: <2016-01-04 14:43:04 t-nissie>
	! Author: Takeshi NISHIMATSU
	!!
	#if defined HAVE_CONFIG_H
	# include "config.h"
	#endif
	program feram_fftw_SoA_out_len
	use, intrinsic :: iso_c_binding
	implicit none
	include 'fftw3.f03'
	integer :: N_TIMES = 100
	integer :: Lx = 64
	integer :: Ly = 64
	integer :: Lz = 64
	integer :: padding_y = 1
	integer :: length_structure = 1
	integer :: flags = FFTW_MEASURE
	character(len=100) :: str = 'FFTW_MEASURE'
	integer :: OMP_GET_MAX_THREADS
	integer :: narg, ireturn
	type(C_PTR) :: plan_r2c_1_out = c_null_ptr
	type(C_PTR) :: plan_c2r_1_out = c_null_ptr

	ireturn = fftw_init_threads()
	call fftw_plan_with_nthreads(OMP_GET_MAX_THREADS())

	ireturn = fftw_import_wisdom_from_filename(C_CHAR_"wisdom"//C_NULL_CHAR)
	if (ireturn.eq.1) then
	write(0,'(a,a,i2,a)') __FILE__, ':', __LINE__, ': Successfully imported FFTW wisdom in current directory.'
	else
	write(0,'(a,a,i2,a)') __FILE__, ':', __LINE__, ': Failed to import FFTW wisdom in current directory.'
	end if

	narg = command_argument_count()
	if (narg.eq.0) then
	! default values
	else if (narg.eq.7) then
	call get_command_argument(1,str); read(str,*) N_TIMES
	call get_command_argument(2,str); read(str,*) Lx
	call get_command_argument(3,str); read(str,*) Ly
	call get_command_argument(4,str); read(str,*) Lz
	call get_command_argument(5,str); read(str,*) padding_y
	call get_command_argument(6,str); read(str,*) length_structure
	call get_command_argument(7,str)
	select case(str)
	case( 'FFTW_ESTIMATE')
	flags = FFTW_ESTIMATE
	case( 'FFTW_PATIENT')
	flags = FFTW_PATIENT
	case( 'FFTW_MEASURE')
	flags = FFTW_MEASURE
	case( 'FFTW_EXHAUSTIVE')
	flags = FFTW_EXHAUSTIVE
	case default
	write(0, '(a,a,i2,a)') __FILE__, ':', __LINE__, ': No such FFTW flag.'
	stop 2
	end select
	else
	write(0, '(a,a,i2,a)') __FILE__, ':', __LINE__, ': Illegal number of arguments.'
	stop 1
	end if
	write(0,'(a,a,i2,a)') __FILE__, ':', __LINE__, ': flags = '//trim(str)

	call feram_fftw_SoA_out(length_structure, N_TIMES, Lx, Ly, Lz, padding_y, flags, &
	& plan_r2c_1_out, &
	& plan_c2r_1_out)

	#ifndef WITH_FFT_MKL
	ireturn = fftw_export_wisdom_to_filename(C_CHAR_"wisdom_SoA"//C_NULL_CHAR)
	#endif
	call fftw_destroy_plan(plan_r2c_1_out)
	call fftw_destroy_plan(plan_c2r_1_out)
	call fftw_cleanup_threads()
	end program feram_fftw_SoA_out_len
	# --Makefile-- for feram_fftw_SoA_out_len
	##
	FC=gfortran
	FFLAGS=-Wall -Ofast -fopenmp -ffree-form
	CPPFLAGS=-I../api
	LDFLAGS=../.libs/libfftw3.a ../threads/.libs/libfftw3_omp.a
	feram_fftw_SoA_out_len: feram_fftw_SoA_out_len.o feram_fftw_SoA_out.o feram_fftw_wisdom_timing.o
	$(FC) $(FFLAGS) -o $@ $^ $(LDFLAGS)
	clean:
	rm -f feram_fftw_*.o feram_fftw_SoA_out_len