MPIを使った並列化の実践(その1)

　某所でMPI(Message Passing Interface)の講習を受けてきた。話を聞いていたら、自分の目的を達成するには、ほんとに基本的な構文だけ知っていればOKのような気がした。

　少なくとも今の俺には、次の問題が解ければいい：


ここで、gは非負の関数とする。さらに、Θはコンパクト。具体例としては、非線形最小二乗法。

　この場合、最上層において、つまりθについてMPIを使って仕事を振り分けることが望ましい（もちろん、g(θ)の計算時間が莫大であることを前提とする）。完全に並列化可能で、ほぼ分散処理となる。
　具体的なコードを次のように書いてみた：

module INCMPI
implicit none
include 'mpif.h'
integer istatus(mpi_status_size)
integer nprocs, myrank, ista, iend, ierror
end
-----------------------------------------------------------------
program MPI
use INCMPI
implicit none
integer, parameter :: N=900
real, dimension(N) :: Theta
real, dimension(:), allocatable :: Objv
real, parameter :: theta_l=0.5, theta_u=1.0
real c, d, obj, int_r
integer i, arg, dummy, int_i

! the space of parameter
Theta(1)=theta_l
Theta(N)=theta_u
c=(theta_u-theta_l)/(N-1)
do i=2,N-1
    Theta(i)=Theta(i-1)+c
enddo

! multi-processing through MPI
call mpi_init(ierror)
call mpi_comm_size(mpi_comm_world,nprocs,ierror)
call mpi_comm_rank(mpi_comm_world,myrank,ierror)

allocate(Objv(nprocs))
int_r=N/nprocs
int_i=floor(int_r)
if (myrank+1<nprocs) then
    ista=int_i*myrank+1
    iend=int_i*(myrank+1)
else
    ista=int_i*myrank+1
    iend=int_i*(myrank+1)+mod(N,nprocs)
endif

! minimization
obj=99999
do i=ista, iend
    call f(Theta(i),d)
    if (d<=min) then
        obj=d
        arg=i
    endif
enddo
Objv(myrank)=obj

! communication
call mpi_allgather(mpi_in_place, dummy, dummy, Objv, 1, mpi_real, mpi_comm_world, ierror)

! summarizing the result
if (obj==min(Objv)) then
    print *, 'The optimal parameter is', Theta(arg),'.'
    print *, 'The minimized value of the objective function is', obj,'.'
endif

deallocate(Objv)
call mpi_finalize(ierror)

end program MPI
-----------------------------------------------------------------
subroutine f(theta,d)
real theta, d
d = g(theta)
end subroutine f
-----------------------------------------------------------------
real, function g(theta)
real theta
g = ... ! nonnegative real number
end function g

ただし、このコードはまだデバッグしてない。まだMPIを使える環境にないため、こんなコードでよかろうに、という想像の世界でございます。

なお、''mpi_in_place''を使っているので、MPI-2に対応している必要がある。