配列代入のLOOP構造分類

Fortranは配列操作が得意なわけですが、配列に要素を計算して代入するときに、その要素が配列のインデックスにどの程度依存するかで並列化に適したLOOP構造を選べます。

配列の要素が全くインデックスによらない場合、インデックスを省略した記法が使えます。(例：b = a)
配列の要素が、その要素のインデックスのみに依るのならば、各要素は独立なので、順序を気にせず並列に計算できます。この場合は、FORALLないしDO CONCURRENTが使えます。FORALLは基本的に代入操作にしか使えず、かつ設計仕様が複雑すぎて実装が難しく、並列化も効かせにくいので、使わないほうがよいとされていますが、１行でLOOP構造が書けるという便利な点もあります。（例：FORALL(i=1:n) a(i) = i**2)DO CONCURRENTはこの欠点を踏まえてFORALLを代替するものとしてFortran2008で導入された命令です。代入操作以外にも使え、並列化にも有利ですが、１行でLOOPを書くことはできません。
最後に配列の要素が、漸化式のように自己のインデックス以外の配列要素にアクセスしている場合は、計算の順序に依存性があるので、基本的にシリアルに計算する必要があります。この場合は従来のDO LOOPを使用することになります。

配列要素に何らかの条件を課す時、前二者ではMASKを利用できます。普通のDO LOOPではIF文を用いて記述することになります。

WHERE (a>0.0) b = LOG(a)
FORALL (i=1:n, a>0.0) b(i) = a(i) * log(real(i)) あるいは DO CONCURRENT (i=1:n, a>0.0); b(i) = a(i) * log(real(i)); END DO
DO i = 1, n ; if (a(i) > 0.0) b(i) = b(i - 1) + LOG(a(i)); END DO

これらとELEMENTALな副プログラムを組み合わせることにより、抽象度が高くコンパイラが並列化しやすいプログラムを書くことが可能になります。

ソースプログラム

program LOOP
  implicit none
  integer, parameter :: n = 5
  real :: a(n), b(n), c(n), d(n)
  integer :: i
  a = [(real(i - (n + 1) / 2),   i = 1, n)]
! no dependence on index
  b = a**2
! local dependence on index   i
  forall(i = 1:n) c(i) = i * a(i)
! local dependence on index   ii
  do concurrent (i = 1:n)
    c(i) = i * a(i)   
  end do
! non-local dependence on index
  d(1) = 0.0
  do i = 2, n
   d(i) = d(i - 1) + a(i)
  end do
  print *, a
  print *
  print *
  print *, b
  print *
  print *, c
  print *
  print *, d
  print *
!
! masking 
!
  b = 0.0
  c = 0.0
  d = 0.0
! no dependence on index
  where(a >= 0.0) b = a**2
! local dependence on index    i
  forall(i = 1:n,  a(i) >= 0.0) c(i) = i * a(i)
! local dependence on index    ii
  do concurrent (i = 1:n,  a(i) >= 0.0)
    c(i) = i * a(i)   
  end do
! non-local dependence on index
  d(1) = 0.0
  do i = 2, n
    if (a(i) >= 0.0) d(i) = d(i - 1) + a(i)
  end do
  print *
  print *, b
  print *
  print *, c
  print *
  print *, d
  print *
!
!  masking with else
!  
  b = 0.0
  c = 0.0
  d = 0.0
! no dependence on index
  where(a >= 0.0) 
    b = a**2
  else where
    b = 0.0 / 0.0 ! NaN
  end where
! local dependence on index   i
  forall(i = 1:n, a(i) >= 0.0) c(i) = i * a(i)
  forall(i = 1:n, a(i) <  0.0) c(i) = 0.0 / 0.0
! local dependence on index   ii
  do concurrent (i = 1:n) 
    if (a(i) >= 0.0) then
      c(i) = i * a(i)
    else
      c(i) = 0.0 / 0.0 ! NaN
    end if  
  end do
! non-local dependence on index
  d(1) = 0.0
  do i = 2, n
    if (a(i) >= 0.0) then
      d(i) = d(i - 1) + a(i)
    else
      d(i) = a(i)
    end if
  end do
  print *
  print *, a
  print *
  print *, b
  print *
  print *, c
  print *
  print *, d
  print *
!  
  stop
end program LOOP

fortran66のブログ

fortran について書きます。

ソースプログラム

実行結果