fortran66のブログ

fortran について書きます。

CAF と MPI 例題比較

W. Gropp 等の Using MPI の最初の例題を CAF で書いて見ます。

つらつら本を読んでいますと、並列プログラミングの場合、1)通信、2)同期、3)排他を考えなければならないようです。通信がなければ単なる独立な実行と同じになります。通信する場合は順序(タイミング)が定まらねばなりませんが、順序が一意に確定していれば同期によって実現します。順序に任意性がある場合はさらに排他処理が必要になります。


プログラム

CAF の場合単方向通信なので、通信と同期が分離していて、明示的に同期を取らなければなりません。MPI の場合基本がブロッキング型の双方向通信なので、通信が暗黙に同期を含んでいます。

また Fortran 2008 の CAF の場合は集団通信用の命令がないため*1、ギャザー・スキャター型の、順序に任意性のある演算の場合、明示的に排他処理をしなければなりません。 MPI は沢山の集団通信用のルーチンを用意して、排他処理を明示的にしなくてよいようにしてあるようです。

CAF
    program CAF3_2
      implicit none
      integer, parameter :: kd = kind(1.0d0)
      real(kd), parameter :: pi_kd = 4 * atan(1.0_kd)
      real(kd) :: pi[*], mypi, h, sum, x, f, a
      integer :: n[*], i, me, ne
      !
      f(a) = 4.0_kd / (1.0_kd + a * a)
      !
      me = this_image()
      ne = num_images()
      do 
        ! Broadcast n  
        if (me == 1) then 
          print *, 'enter the number of intervals: (0 quits) '
          read *, n
          pi = 0.0_kd
          sync images(*)
        else
          sync images(1)
          n = n[1]
        end if
        if (n == 0) exit
        !
        h = 1.0_kd / n
        sum = 0.0_kd
        do i = me, n, ne
          x = h * (real(i, kd) - 0.5_kd)
          sum = sum + f(x)
        end do
        mypi = h * sum
        !
        ! reduction
        critical 
          pi[1] = pi[1] + mypi
        end critical
        sync all
        if (me == 1) then
          print *, 'pi is', pi, ' Error is', abs(pi - pi_kd)
        end if  
      end do
    end program CAF3_2
CAF 実行結果
 enter the number of intervals: (0 quits)
10
 pi is   3.14242598500110       Error is  8.333314113051493E-004
 enter the number of intervals: (0 quits)
100
 pi is   3.14160098692312       Error is  8.333333331833614E-006
 enter the number of intervals: (0 quits)
1000
 pi is   3.14159273692313       Error is  8.333333401111531E-008
 enter the number of intervals: (0 quits)
10000
 pi is   3.14159265442313       Error is  8.333342904620622E-010
 enter the number of intervals: (0 quits)
0
MPI
    program mpi3_2
      use mpi
      implicit none
      integer, parameter :: kd = kind(0.0d0)
      real(kd), parameter :: pi25dt = 3.141592653589793238462643_kd
      real(kd) :: mypi, pi, h, sum, x, f, a
      integer :: i, n, ierr, myid, numprocs
      ! statement function
      f(a) = 4.0_kd / (1.0_kd + a * a)
      !
      call MPI_INIT(ierr)
      call MPI_COMM_RANK(MPI_COMM_WORLD, myid, ierr)
      call MPI_COMM_SIZE(MPI_COMM_WORLD, numprocs, ierr)      
      do 
        if (myid == 0) then
          print *, 'Enter the number of intervals: (0 quits) '
          read  *, n
        end if    
      ! broadcast n
        call MPI_BCAST(n, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, ierr)
      !
        if (n == 0) exit
      !
        h = 1.0_kd / n
        sum = 0.0_kd
        do i = myid + 1, n, numprocs
          x = h * (real(i, kd) - 0.5_kd) 
          sum = sum + f(x)
        end do  
        mypi = h * sum
      !  collect all the partial sums  
        call MPI_REDUCE(mypi, pi, 1, MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, ierr)
      !
        if (myid == 0) then
          print *, 'pi is ', pi, ' Error is', abs(pi - pi25dt)
        end if
      end do  
      call MPI_FINALIZE(ierr)
    end program mpi3_2
MPI 実行結果
C:\Users\O\Documents\Visual Studio 2015\Projects\MPI\MPI\x64\Debug>mpiexec -n 4 mpi.exe
 Enter the number of intervals: (0 quits)
10
 pi is    3.14242598500110       Error is  8.333314113051493E-004
 Enter the number of intervals: (0 quits)
100
 pi is    3.14160098692312       Error is  8.333333331833614E-006
 Enter the number of intervals: (0 quits)
1000
 pi is    3.14159273692313       Error is  8.333333312293689E-008
 Enter the number of intervals: (0 quits)
10000
 pi is    3.14159265442312       Error is  8.333307377483834E-010
 Enter the number of intervals: (0 quits)
0

*1:Fortran 2015 で導入予定です。