t-nissieの日記： 【電脳】PS3 (Cell) ってすごかったんだなぁ
日記 by
t-nissie
PS3 (Cell) ってすごかったんだなぁ（過去形）。
よいコンパイラがあってプログラミングが簡単だったなら生き残れたかなぁ。
単精度３次元FFT対決 cuFFT (GPU) vs FFTW (x86_64) vs FFTW (Cell)
t-nissieさんのトモダチの日記、みんなの日記も見てね。 スラッシュドットのストーリーを選ぶための補助をお願いします。
t-nissieの日記： 【電脳】円高だしFree Software Foundationに寄付をしよう！
日記 by
t-nissie
円高だしFree Software Foundationに寄付をしよう！
クレジットカードも使えるよ。
t-nissieの日記： 【電脳】GPU/CUDAプログラミング cuFFTを試したい 1
日記 by
t-nissie
５年以上前に買った小型のパソコンHP dc5100 (Pentium 4) に
玄人指向のGF-GT520-LE1GHを入れればUbuntuで
GPU/CUDAプログラミングが体験できるようになるのだろうか。
とりあえずcuFFTを試したい。
HP dc5100 http://h50146.www5.hp.com/products/desktops/old/dc5100sf_ct/pentium4_model.html
GF-GT520-LE1GH http://www.kuroutoshikou.com/modules/display/?iid=1598
t-nissieの日記： 【電脳】FFTWのOpenMP並列でいまいち３次元r2c, c2rの速さが出ない【その３】
日記 by
t-nissie
$ for i in `jot 3`; do OMP_NUM_THREADS=8 ./r2c_3d_test_omp 32 32 160 100000 8; done ; for i in `jot 3`; do OMP_NUM_THREADS=8 ./r2c_3d_test_threads 32 32 160 100000 8; done
Lx = 32, Ly = 32, Lz = 160, N = 163840, M = 100000, NTHREADS = 8
Address of r = 0x0000000100400000
Address of c = 0x0000000100540000
FFT starts
FFT ends 256.37600
Lx = 32, Ly = 32, Lz = 160, N = 163840, M = 100000, NTHREADS = 8
Address of r = 0x0000000100400000
Address of c = 0x0000000100540000
FFT starts
FFT ends 254.20400
Lx = 32, Ly = 32, Lz = 160, N = 163840, M = 100000, NTHREADS = 8
Address of r = 0x0000000100400000
Address of c = 0x0000000100540000
FFT starts
FFT ends 263.79200
Lx = 32, Ly = 32, Lz = 160, N = 163840, M = 100000, NTHREADS = 8
Address of r = 0x0000000100400000
Address of c = 0x0000000100540000
FFT starts
FFT ends 239.95000
Lx = 32, Ly = 32, Lz = 160, N = 163840, M = 100000, NTHREADS = 8
Address of r = 0x0000000100400000
Address of c = 0x0000000100540000
FFT starts
FFT ends 237.44200
Lx = 32, Ly = 32, Lz = 160, N = 163840, M = 100000, NTHREADS = 8
Address of r = 0x0000000100400000
Address of c = 0x0000000100540000
FFT starts
FFT ends 236.12100
$
! r2c_3d_test.F -*-f90-*-
! Time-stamp: <2011-11-29 15:53:27 takeshi>
! Author: Takeshi NISHIMATSU
!!
#if defined(__PGI) || defined(SR11000) || defined(__sparc)
# define command_argument_count iargc
# define get_command_argument getarg
#endif
program r2c_3d_test
implicit none
real*8, allocatable :: r(:,:,:)
complex*16, allocatable :: c(:,:,:)
integer*8 :: plan_r2c, plan_c2r, address, count0, count1, count_rate
character(len=30) :: str
integer :: Lx, Ly, Lz, M, NTHREADS, i, j, ireturn
real*8 :: N_inv
# include "fftw3.f"
call get_command_argument(1,str); read(str,*) Lx
call get_command_argument(2,str); read(str,*) Ly
call get_command_argument(3,str); read(str,*) Lz
call get_command_argument(4,str); read(str,*) M
call get_command_argument(5,str); read(str,*) NTHREADS
write(6,'(3(a,i4),2(a,i7),a,i3)') 'Lx = ', Lx, &
& ', Ly = ', Ly, &
& ', Lz = ', Lz, &
& ', N = ', Lx*Ly*Lz, &
& ', M = ', M, &
& ', NTHREADS = ', NTHREADS
N_inv = 1.0d0 / Lx / Ly / Lz
call dfftw_init_threads(ireturn)
call dfftw_plan_with_nthreads(NTHREADS)
allocate(r(0:Lx-1, 0:Ly-1, 0:Lz-1))
allocate(c(0:Lx/2, 0:Ly-1, 0:Lz-1))
write (*,'(a,z16.16)') 'Address of r = 0x', address(r) ! Check 16-bit alignment,
write (*,'(a,z16.16)') 'Address of c = 0x', address(c) ! or SSE2 won't be used.
call dfftw_plan_dft_r2c_3d(plan_r2c, Lx, Ly, Lz, r, c, FFTW_PATIENT)
call dfftw_plan_dft_c2r_3d(plan_c2r, Lx, Ly, Lz, c, r, FFTW_PATIENT)
r(:,:,:) = 0.1d0
write(6,'(a)') 'FFT starts'
call flush(6)
call system_clock(count0)
do i = 1, M
call dfftw_execute(plan_r2c)
!$omp parallel do
do j = 0, Lz-1
c(:,:,j) = c(:,:,j) * N_inv
end do
!$omp end parallel do
call dfftw_execute(plan_c2r)
end do
call system_clock(count1, count_rate)
write(6,'(a,f10.5)') 'FFT ends', dble(count1-count0)/count_rate
call dfftw_cleanup_threads(ireturn)
end program r2c_3d_test
!Local variables:
! compile-command: "make -k && ./r2c_3d_test_omp 32 32 160 100 4"
!End:
t-nissieの日記： 【電脳】FFTWのOpenMP並列でいまいち３次元r2c, c2rの速さが出ない【その２】 1
日記 by
t-nissie
以下のようなプログラムでいろいろ試しています。
topで見るかぎりあまり並列化されていないようです。
フーリエ変換する三次元配列のサイズが小さいせいか
あまり並列化しないほうが速いのかもしれません。
FFTW_PATIENTだとdfftw_plan_dft_r2c_3dなどで
planを作るのに数分かかります。以下で試したところ
fftw3_threadsライブラリを使うよりfftw3_ompを
使うほうが速いようです。16-bit alignmentはされて
いるようなのでSSE2は使われているものと思います。
どんなplanになったかを表示する関数はないのかしらん。
fftw_describe(plan), fftw_show(plan)みたいな。
これから、OpenMPでネストして速くなるか検証する
予定。
$ for i in `jot 9`; do OMP_NUM_THREADS=4 ./r2c_3d_test_omp 32 32 160 1000 4 | grep ends; done ; \
for i in `jot 9`; do OMP_NUM_THREADS=4 ./r2c_3d_test_threads 32 32 160 1000 4 | grep ends; done
FFT ends 3.17200
FFT ends 3.02300
FFT ends 3.19100
FFT ends 3.06800
FFT ends 3.23800
FFT ends 3.12100
FFT ends 3.20900
FFT ends 3.24900
FFT ends 3.19000
FFT ends 4.69000
FFT ends 4.56500
FFT ends 4.61900
FFT ends 4.64000
FFT ends 4.63200
FFT ends 4.62900
FFT ends 4.65000
FFT ends 4.60900
FFT ends 4.69600
# -*-Makefile-*- for testing FFTW3
##
FC = gfortran
FFLAGS = -Wall -ffree-form -O3 -pipe -fopenmp -fopenmp
all: r2c_3d_test_threads r2c_3d_test_omp
r2c_3d_test_threads: r2c_3d_test.o address.o
gfortran $(FFLAGS) -lfftw3 -lfftw3_threads -o $@ $^
r2c_3d_test_omp: r2c_3d_test.o address.o
gfortran $(FFLAGS) -lfftw3 -lfftw3_omp -o $@ $^
clean:
rm -f *.o
! r2c_3d_test.F -*-f90-*-
! Time-stamp: <2011-11-28 20:21:46 t-nissie>
!!
#if defined(__PGI) || defined(SR11000) || defined(__sparc)
# define command_argument_count iargc
# define get_command_argument getarg
#endif
program r2c_3d_test
implicit none
real*8, allocatable :: r(:,:,:)
complex*16, allocatable :: c(:,:,:)
integer*8 :: plan_r2c, plan_c2r, address, count0, count1, count_rate
character(len=30) :: str
integer :: Lx, Ly, Lz, M, NTHREADS, i, j, ireturn
real*8 :: N_inv
# include "fftw3.f"
call get_command_argument(1,str); read(str,*) Lx
call get_command_argument(2,str); read(str,*) Ly
call get_command_argument(3,str); read(str,*) Lz
call get_command_argument(4,str); read(str,*) M
call get_command_argument(4,str); read(str,*) NTHREADS
N_inv = 1.0d0 / Lx / Ly / Lz
call dfftw_init_threads(ireturn)
call dfftw_plan_with_nthreads(NTHREADS)
allocate(r(0:Lx-1, 0:Ly-1, 0:Lz-1))
allocate(c(0:Lx/2, 0:Ly-1, 0:Lz-1))
write (*,'(a,z16.16)') 'Address of r = 0x', address(r) ! Check 16-bit alignment,
write (*,'(a,z16.16)') 'Address of c = 0x', address(c) ! or SSE2 won't be used.
call dfftw_plan_dft_r2c_3d(plan_r2c, Lx, Ly, Lz, r, c, FFTW_PATIENT)
call dfftw_plan_dft_c2r_3d(plan_c2r, Lx, Ly, Lz, c, r, FFTW_PATIENT)
r(:,:,:) = 0.1d0
write(6,'(a)') 'FFT starts'
call flush(6)
call system_clock(count0)
do i = 1, M
call dfftw_execute(plan_r2c)
!$omp parallel do
do j = 0, Lz-1
c(:,:,j) = c(:,:,j) * N_inv
end do
!$omp end parallel do
call dfftw_execute(plan_c2r)
end do
call system_clock(count1, count_rate)
write(6,'(a,f10.5)') 'FFT ends', dble(count1-count0)/count_rate
call dfftw_cleanup_threads(ireturn)
end program r2c_3d_test
!Local variables:
! compile-command: "make -k && ./r2c_3d_test 32 32 160 100 4"
!End:
unsigned long address_(void *x)
{
return (unsigned long)x;
}
t-nissieの日記： 【電脳】FFTWのOpenMP並列でいまいち３次元r2c, c2rの速さが出ない 1
日記 by
t-nissie
FFTWのOpenMP並列でいまいち３次元r2c, c2rの速さが出ない。
planを作るとき並列化効率が上がらないと
thread数を指定より勝手に下げてしまうみたい。
そのときネストして並列化してやる必要があるのだろうか？？？
http://www.fftw.org/doc/How-Many-Threads-to-Use_003f.html
t-nissieの日記： 【電脳】Mac OS X Lion (10.7.x) の /usr/bin/xattr って速くなったのかな
日記 by
t-nissie
なんか10.5.xのはものすごく遅かったのだけど。
ってかまだ10.5.x使ってる。
t-nissieの日記： 【野球】１毛差でクライマックスシリーズを逃すオリックスって…
日記 by
t-nissie
１毛差でクライマックスシリーズを逃すオリックスって…
一方、田口壮外野手は肩を手術 http://taguchiso.blogspot.com/
t-nissieの日記： 【電脳】どっちが速いのか 7
日記 by
t-nissie
こういうのって、どっちが速いのかしらん。
Nxの大きさとキャッシュの大きさに関係するのかしらん。
プログラム全体からしたらゴミみたいなところだから影響は小さいのだけど。
equivalence文を使えばもっとクールに書けるかな。
before:
subroutine dfftw_plan_dft_c2r_3d(plan, Nx, Ny, Nz, c, r, flags)
implicit none
integer*8, intent(in) :: plan
integer, intent(in) :: Nx, Ny, Nz
complex*16 :: c(0:Nx/2, 0:Ny-1, 0:Nz-1) !CAUTION "/2"!
real*8, intent(out):: r(0:Nx-1, 0:Ny-1, 0:Nz-1)
integer, intent(in) :: flags
integer ix, iy, iz, icon
call DM_V3DRCF(c, 2*(Nx/2+1), Nx, Ny, Nz, 1, -1, ICON)
do iz = 0, Nz-1
do iy = 0, Ny-1
do ix = 0, Nx-1
if (mod(ix,2)==0) then
r(ix,iy,iz) = dble(c(ix/2,iy,iz))
else
r(ix,iy,iz) = aimag(c(ix/2,iy,iz))
end if
end do
end do
end do
end subroutine dfftw_plan_dft_c2r_3d
after:
subroutine dfftw_plan_dft_c2r_3d(plan, Nx, Ny, Nz, c, r, flags)
implicit none
integer*8, intent(in) :: plan
integer, intent(in) :: Nx, Ny, Nz
complex*16 :: c(0:Nx/2, 0:Ny-1, 0:Nz-1) !CAUTION "/2"!
real*8, intent(out):: r(0:Nx-1, 0:Ny-1, 0:Nz-1)
integer, intent(in) :: flags
integer ix, iy, iz, icon
call DM_V3DRCF(c, 2*(Nx/2+1), Nx, Ny, Nz, 1, -1, ICON)
do iz = 0, Nz-1
do iy = 0, Ny-1
do ix = 0, Nx-1, 2
r(ix,iy,iz) = dble(c(ix/2,iy,iz))
end do
do ix = 1, Nx-1, 2
r(ix,iy,iz) = aimag(c(ix/2,iy,iz))
end do
end do
end do
end subroutine dfftw_plan_dft_c2r_3d
t-nissieの日記： 【電脳】ログインしている状態だとコメントの数が表示されない 3
日記 by
t-nissie
トップページで