t-nissieの日記： 【電脳】どっちが速いのか 7
こういうのって、どっちが速いのかしらん。
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
そもそも ix に対する if 文がループの一番内側って… (スコア:1)
それを一番外に出しておけば、「1つの条件判定内で、ループ」になるから、ずっと最適化しやすくなると思うんだが?
fjの教祖様
Re:そもそも ix に対する if 文がループの一番内側って… (スコア:1)
やりたいことは，配列をごそっとコピーすることです．
ただし，Nxが奇数のときに端が余るので，そこを気にしないといけません．
CとFortranとでは配列のメモリへの格納の順番みたいなモノがちがいます．
Fortranでは，ここのixのほうがメモリ上で横に並びます．
それでも条件分岐をループの外にした方が得なのかなぁ，ってのが疑問です．
love && peace && free_software
t-nissie
Re:そもそも ix に対する if 文がループの一番内側って… (スコア:1)
たぶん最大のポイントは、これがコンパイラに伝わっていない、という点でしょう。
その理由は、奇数と偶数で、呼び出す関数が違うという形で書かれていること。そして(これは予想ですが)、おそらくその関数の実態がどのような存在なのか、この部分をコンパイルしているときにはノーヒントだ、というのが問題。
関数コールのたびに Cache やメモリアクセスに関する予測は「予測不能」になってしまう。ループの最も内側で if 文を用意して、それぞれの分岐結果が「予測不能」なブロックで構成されていると、最適化予測はほとんど不可能になります。結果、安全サイドに倒さざるを得なくなる。
最も内側のループを同一性質でできたループ2つに分けると速度が上がるなら、ix によるループをいちばん外側に出すとさらに速度は上がると思います。
でも、それ以前に、ループの一番内側にいる関数2つが、そもそも何なのかをコンパイラに教える方法がないか、考えたほうがいいのは事実です。両方の関数をインライン展開できるだけのヒントがあれば(たとえコンパイル時にインライン展開していないとしても)、それなりの「配慮」ってやつをしてくれるコンパイラは多いです。
fjの教祖様
afterの方が速い (スコア:0)
キャッシュ容量云々よりも、少し前の分岐予測なら一番嫌いな分岐先が毎回変わるパターンになっているからです。
＃１回前の履歴だけしか保持していない場合、毎回ミスするパターンです
in-order処理で分岐予測が無い場合でも、mod命令分だけ速くなるはずです。
Re:afterの方が速い (スコア:1)
このくらいの分岐とmod命令ならコンパイラが空気を読んで最適化してくれないもんですかね．
だめですよね．
love && peace && free_software
t-nissie
#2035419ですが、次が最速じゃね？ (スコア:0)
r(ix,iy,iz) = dble(c(ix/2,iy,iz))
r(ix+1,iy,iz) = aimag(c(ix/2,iy,iz))
end do
if ( ix < Nx ) then
r(ix,iy,iz) = dble(c(ix/2,iy,iz))
end if
Nxは偶数って制約ありそうなので、最後のif文は使われないと思いますが念のため。
(ix判定が不安なら↓みたくNx判定でも可。但しNx==0でエラーですがｗ)
if (mod(Nx,2)==0) then
r(Nx-1,iy,iz) = dble(c((Nx-1)/2,iy,iz))
end if
Re:#2035419ですが、次が最速じゃね？ (スコア:1)
Nxが偶数ならcとrとは同じサイズになるのでずっと単純に書けます．
まあ，プログラム全体としてはNxを奇数にする意味も必要もないのですが…
love && peace && free_software
t-nissie