課題S1解説【Fortran言語編】

RIKEN AICS HPC Summer School 2014
中島研吾（東大・情報基盤センター）
横川三津夫（神戸大・計算科学教育センター）

MPI Programming

課題S1 （1/2）• 「<$P-S1>/a1.0～a1.3」， 「<$P-S1>/a2.0~a2.3」から局
所ベクトル情報を読み込み，全体ベクトルのノルム（||x||）を求めるプログラムを作成する（S1-1）．
– ノルム ||x|| は，各要素の２乗の和の平方根である．
– <$P-S1>file.f，<$P-S1>file2.fをそれぞれ参考にする．

• 「<$P-S1>/a2.0~a2.3」から局所ベクトル情報を読み込み，「全体ベクトル」情報を各プロセッサに生成するプログラムを作成する．MPI_Allgathervを使用する（S1-2）．

2

MPI Programming

課題S1 （2/2）• 下記の数値積分を台形公式によって求めるプログラムを作成する．MPI_Reduce，MPI_Bcast等を使用して並列化を実施
し，プロセッサ数を変化させた場合の計算時間を測定する（S1-3）．
4
0 1  x 2 dx
1

N


1
x  f1  f N 1 
2 f i 
2 
i 2




3

S1-ref

ファイルコピーFORTRANユーザー>$ cd <$P-TOP>>$ cp /tmp/2014summer/F/s1r-f.tar .
>$ tar xvf s1r-f.tar
Cユーザー>$ cd <$P-TOP>>$ cp /tmp/2014summer/C/s1r-c.tar .
>$ tar xvf s1r-c.tar
ディレクトリ確認>$ ls
mpi>$ cd mpi/S1-ref

このディレクトリを本講義では <$P-S1r> と呼ぶ。
<$P-S1r> = <$P-TOP>/mpi/S1-ref4

S1-ref

S1-1：局所ベクトル読み込み，ノルム計算• 「<$P-S1>/a1.0～a1.3」から局所ベクトル情報を読み込み，全体ベクトルのノルム（||x||）を求めるプログラムを作成する（S1-1）。
• MPI_Allreduce（またはMPI_Reduce）の利用• ワンポイントアドバイス– 変数の中身を逐一確認しよう !

S1-1

5

S1-ref

MPI_REDUCE•

P#0 A0 B0 C0 D0P#1 A1 B1 C1 D1P#0 op.A0-A3 op.B0-B3 op.C0-C3 op.D0-D3
Reduce

P#1

P#2 A2 B2 C2 D2

P#2

P#3 A3 B3 C3 D3

P#3

コミュニケーター 「comm」内の，各プロセスの送信バッファ「sendbuf」について，
演算「op」を実施し，その結果を1つの受信プロセス「root」の受信バッファ「recbuf」に格納する。
– 総和，積，最大，最小 他• call MPI_REDUCE(sendbuf,recvbuf,count,datatype,op,root,comm,ierr)
任意
任意

I
O

– count
整数– datatype 整数

II– sendbuf– recvbuf

FORTRANC– op

送信バッファの先頭アドレス，
受信バッファの先頭アドレス，
タイプは「datatype」により決定
メッセージのサイズ
メッセージのデータタイプ

MPI_INTEGER, MPI_REAL, MPI_DOUBLE_PRECISION, MPI_CHARACTER etc.
MPI_INT, MPI_FLOAT, MPI_DOUBLE, MPI_CHAR etc

整数

I

計算の種類

MPI_MAX, MPI_MIN, MPI_SUM, MPI_PROD, MPI_LAND, MPI_BAND etc

ユーザーによる定義も可能： MPI_OP_CREATE– root
整数I
受信元プロセスのID（ランク）– comm
整数
I
コミュニケータを指定する– ierr
整数
O
完了コード

6

S1-ref

送信バッファと受信バッファ• MPIでは「送信バッファ」，「受信バッファ」という変数がしばしば登場する。
• 送信バッファと受信バッファは必ずしも異なった名称の配列である必要はないが，必ずアドレスが異なっていなければならない。

7

S1-ref

MPI_Reduce/Allreduceの “op”
call MPI_REDUCE(sendbuf,recvbuf,count,datatype,op,root,comm,ierr)• MPI_MAX，MPI_MIN• MPI_SUM，MPI_PROD• MPI_LAND

最大値，最小値
総和，積
論理AND

real(kind=8):: X0, X1
call MPI_REDUCE(X0, X1, 1, MPI_DOUBLE_PRECISION, MPI_MAX, 0, <comm>, ierr)real(kind=8):: X0(4), XMAX(4)call MPI_REDUCE(X0, XMAX, 4, MPI_DOUBLE_PRECISION, MPI_MAX, 0, <comm>, ierr)
8

S1-ref

MPI_BCAST•P#0 A0 B0 C0 D0
P#1

P#0 A0 B0 C0 D0BroadcastP#1 A0 B0 C0 D0

P#2

P#2 A0 B0 C0 D0

P#3

P#3 A0 B0 C0 D0

コミュニケーター 「comm」内の一つの送信元プロセス「root」のバッファ「buffer」
から，その他全てのプロセスのバッファ「buffer」にメッセージを送信。

• call MPI_BCAST (buffer,count,datatype,root,comm,ierr)– buffer

任意I/O– count

整数I

バッファの先頭アドレス，
タイプは「datatype」により決定
メッセージのサイズ– datatype 整数

I

メッセージのデータタイプ

FORTRANC– root– comm– ierr

MPI_INTEGER, MPI_REAL, MPI_DOUBLE_PRECISION, MPI_CHARACTER etc.
MPI_INT, MPI_FLOAT, MPI_DOUBLE, MPI_CHAR etc.

整数
整数
整数

I
I
O

送信元プロセスのID（ランク）
コミュニケータを指定する
完了コード

9

S1-ref

MPI_ALLREDUCE

P#0 op.A0-A3 op.B0-B3 op.C0-C3 op.D0-D3

P#0 A0 B0 C0 D0P#1 A1 B1 C1 D1All reduceP#1 op.A0-A3 op.B0-B3 op.C0-C3 op.D0-D3

P#2 A2 B2 C2 D2P#2 op.A0-A3 op.B0-B3 op.C0-C3 op.D0-D3

P#3 A3 B3 C3 D3P#3 op.A0-A3 op.B0-B3 op.C0-C3 op.D0-D3• MPI_REDUCE + MPI_BCAST• 総和，最大値を計算したら，各プロセスで利用したい場合が多い• call MPI_ALLREDUCE(sendbuf,recvbuf,count,datatype,op, comm,ierr)– sendbuf– recvbuf

任意
任意

I
O

– count– datatype– op– comm– ierr

整数
整数
整数
整数
整数

I
I
I
I
O

送信バッファの先頭アドレス，
受信バッファの先頭アドレス，
タイプは「datatype」により決定
メッセージのサイズ
メッセージのデータタイプ
計算の種類
コミュニケータを指定する
完了コード
10

S1-ref

S1-1：局所ベクトル読み込み，ノルム計算
均一長さベクトルの場合（a1.*）: s1-1-for_a1.f/c
implicit REAL*8 (A-H,O-Z)
include 'mpif.h'
integer :: PETOT, my_rank, SOLVER_COMM, ierrreal(kind=8), dimension(8) :: VEC
character(len=80):: filename
call MPI_INIT(ierr)
call MPI_COMM_SIZE (MPI_COMM_WORLD, PETOT, ierr )call MPI_COMM_RANK (MPI_COMM_WORLD, my_rank, ierr )if (my_rank.eq.0) filename= 'a1.0'
if (my_rank.eq.1) filename= 'a1.1'
if (my_rank.eq.2) filename= 'a1.2'
if (my_rank.eq.3) filename= 'a1.3'

write(filename,'(a,i1.1)') 'a1.', my_rank

N=8open (21, file= filename, status= 'unknown')
do i= 1, N
read (21,*) VEC(i)enddosum0= 0.d0
do i= 1, N
sum0= sum0 + VEC(i)**2
enddo
call MPI_allREDUCE (sum0, sum,
sum= dsqrt(sum)

中身を書き出して見よう

call MPI_Allreduce(sendbuf,recvbuf,count,datatype,op, comm,ierr)
1, MPI_DOUBLE_PRECISION, MPI_SUM, MPI_COMM_WORLD, ierr)

if (my_rank.eq.0) write (*,'(1pe16.6)') sum
call MPI_FINALIZE (ierr)
stop
end

S1-1

11

S1-ref

S1-1：局所ベクトル読み込み，ノルム計算
不均一長さベクトルの場合（a2.*）：s1-1-for_a2.f/c
implicit REAL*8 (A-H,O-Z)
include 'mpif.h'
integer :: PETOT, my_rank, SOLVER_COMM, ierrreal(kind=8), dimension(:), allocatable :: VEC, VEC2
character(len=80):: filename
call MPI_INIT(ierr)
call MPI_COMM_SIZE (MPI_COMM_WORLD, PETOT, ierr )call MPI_COMM_RANK (MPI_COMM_WORLD, my_rank, ierr )if (my_rank.eq.0) filename= 'a2.0'
if (my_rank.eq.1) filename= 'a2.1'
if (my_rank.eq.2) filename= 'a2.2'
if (my_rank.eq.3) filename= 'a2.3'
open (21, file= filename, status= 'unknown')
read (21,*) N
allocate (VEC(N))
do i= 1, N
read (21,*) VEC(i)enddosum0= 0.d0
do i= 1, N
sum0= sum0 + VEC(i)**2
enddo
call MPI_allREDUCE (sum0, sum,
sum= dsqrt(sum)

中身を書き出して見よう

call MPI_Allreduce(sendbuf,recvbuf,count,datatype,op, comm,ierr)
1, MPI_DOUBLE_PRECISION, MPI_SUM, MPI_COMM_WORLD, ierr)

if (my_rank.eq.0) write (*,'(1pe16.6)') sum

S1-1

call MPI_FINALIZE (ierr)
stop
end
12

S1-ref

実 行（課題S1-1）
FORTRAN$ cd <$P-S1r>$ mpifrtpx –Kfast s1-1-for_a1.f$ mpifrtpx –Kfast s1-1-for_a2.f(modify “go4.sh”)$ pjsub go4.shC$ cd <$P-S1r>$ mpifccpx –Kfast s1-1-for_a1.c$ mpifccpx –Kfast s1-1-for_a2.c(modify “go4.sh”)$ pjsub go4.sh
13

S1-ref

S1-1：局所ベクトル読み込み，ノルム計算
計算結果
予め求めておいた答え
a1.* 1.620882475690326E+03
a2.* 1.222184928723964E+03$> f95 dot-a1.f$> ./a.out$> f95 dot-a2.f$> ./a.out

計算結果
a1.* 1.620882475690326E+03
a2.* 1.222184928723964E+03
S1-1

14

S1-ref

S1-1：局所ベクトル読み込み，ノルム計算
SENDBUFとRECVBUFを同じにしたら・・・
正
call MPI_allREDUCE(sum0, sum, 1, MPI_DOUBLE_PRECISION,MPI_SUM, MPI_COMM_WORLD, ierr)誤
call MPI_allREDUCE(sum0, sum0, 1, MPI_DOUBLE_PRECISION,MPI_SUM, MPI_COMM_WORLD, ierr)

S1-1

15

S1-ref

S1-1：局所ベクトル読み込み，ノルム計算
SENDBUFとRECVBUFを同じにしたら・・・
正
call MPI_allREDUCE(sum0, sum, 1, MPI_DOUBLE_PRECISION,MPI_SUM, MPI_COMM_WORLD, ierr)誤
call MPI_allREDUCE(sum0, sum0, 1, MPI_DOUBLE_PRECISION,MPI_SUM, MPI_COMM_WORLD, ierr)正
call MPI_allREDUCE(sumK(1), sumK(2), 1, MPI_DOUBLE_PRECISION,MPI_SUM, MPI_COMM_WORLD, ierr)

これバッファが重なっていないのでOK

S1-1

16

S1-ref

S1-2：局所ベクトルから全体ベクトル生成• 「<$P-S1>/a2.0~a2.3」から局所ベクトル情報を読み込み，「全体ベクトル」情報を各プロセッサに生成するプログラムを作成する。MPI_Allgathervを使用する。

S1-2

17

S1-ref

S1-2：局所ベクトルから全体ベクトル生成
ALLGATHERVを使う場合（1/5）

PE#0

PE#0

PE#1

PE#1

PE#2

PE#2
MPI_Allgatherv

PE#3

S1-2

PE#3

18

S1-ref

MPI_ALLGATHERV• MPI_ALLGATHER の可変長さベクトル版– 「局所データ」から「全体データ」を生成する• call MPI_ALLGATHERV (sendbuf, scount, sendtype, recvbuf,
rcounts, displs, recvtype, comm, ierr)
–
–
–
–
–
–
–
–
–

S1-2

sendbuf 任意
scount整数
sendtype 整数
recvbuf 任意
rcounts 整数
displs整数
recvtype 整数
comm
整数
ierr
整数

I
I
I
O
I
I
I
I
O

送信バッファの先頭アドレス，
送信メッセージのサイズ
送信メッセージのデータタイプ
受信バッファの先頭アドレス，
受信メッセージのサイズ（配列：サイズ＝PETOT）
受信メッセージのインデックス（配列：サイズ＝PETOT+1）
受信メッセージのデータタイプ
コミュニケータを指定する
完了コード

19

S1-ref

MPI_ALLGATHERV（続き）• call MPI_ALLGATHERV (sendbuf, scount, sendtype, recvbuf,
rcounts, displs, recvtype, comm, ierr)
整数 I 受信メッセージのサイズ（配列：サイズ＝PETOT）
整数 I 受信メッセージのインデックス（配列：サイズ＝PETOT+1）– rcounts– displs– この2つの配列は，最終的に生成される「全体データ」のサイズに関する配列であるため，各プロセスで配列の全ての値が必要になる：• もちろん各プロセスで共通の値を持つ必要がある。

– 通常はstride(i)=rcounts(i)

PE#0

PE#1

stride(1)rcounts(1)displs(1)=0

S1-2

stride(2)rcounts(2)PE#2

PE#(m-2) PE#(m-1)

stride(3)stride(m-1)

stride(m)rcounts(3)rcounts(m-1)

rcounts(m)displs(2)=displs(1) + stride(1)size(recvbuf)= displs(PETOT+1)= sum(stride)

displs(m+1)=
displs(m) + stride(m)20S1-ref

PE#0

PE#2

N

N

rcounts(3)

PE#1

N

stride(1)displs(2)rcounts(2)局所データから全体データを生成する

rcounts(1)MPI_ALLGATHERV
でやっていること

displs(1)stride(2)displs(3)stride(3)displs(4)Nrcounts(4)

PE#3

stride(4)displs(5)S1-2

局所データ：sendbuf

全体データ：recvbuf

21

S1-ref

MPI_ALLGATHERV
でやっていること
PE#0

PE#1

N

N

displs(2)stride(2)= rcounts(2)displs(3)stride(3)= rcounts(3)rcounts
(4)

PE#3

N

stride(1)= rcounts(1)rcounts(3)

PE#2

N

rcounts(1) rcounts(2)局所データから全体データを生成する

displs(1)stride(4)= rcounts(4)displs(4)displs(5)S1-2

局所データ：sendbuf

全体データ：recvbuf

22

S1-ref

MPI_ALLGATHERV詳細（1/2）• call MPI_ALLGATHERV (sendbuf, scount, sendtype, recvbuf,
rcounts, displs, recvtype, comm, ierr)– rcounts– displs

整数
整数

I
I

受信メッセージのサイズ（配列：サイズ＝PETOT）
受信メッセージのインデックス（配列：サイズ＝PETOT+1）• rcounts– 各PEにおけるメッセージサイズ：局所データのサイズ• displs– 各局所データの全体データにおけるインデックス– displs(PETOT+1)が全体データのサイズ
PE#0

PE#1

stride(1)rcounts(1)displs(1)=0

S1-2

stride(2)rcounts(2)PE#2

PE#(m-2) PE#(m-1)

stride(3)stride(m-1)

stride(m)rcounts(3)rcounts(m-1)

rcounts(m)displs(2)=displs(1) + stride(1)size(recvbuf)= displs(PETOT+1)= sum(stride)

displs(m+1)=
displs(m) + stride(m)23S1-ref

MPI_ALLGATHERV詳細（2/2）• rcountsとdisplsは各プロセスで共通の値が必要– 各プロセスのベクトルの大きさ N をallgatherして，rcountsに相当するベクトルを作る。
– rcountsから各プロセスにおいてdisplsを作る（同じものができる）。
• stride(i)= rcounts(i) とする– rcountsの和にしたがってrecvbufの記憶領域を確保する。
PE#0

PE#1

stride(1)rcounts(1)displs(1)=0

S1-2

stride(2)rcounts(2)PE#2

PE#(m-2) PE#(m-1)

stride(3)stride(m-1)

stride(m)rcounts(3)rcounts(m-1)

rcounts(m)displs(2)=displs(1) + stride(1)size(recvbuf)= displs(PETOT+1)= sum(stride)

displs(m+1)=
displs(m) + stride(m)24S1-ref

MPI_ALLGATHERV使用準備
例題：<$P-S1>/agv.fc• “a2.0”~”a2.3”から，全体ベクトルを生成する。
• 各ファイルのベクトルのサイズが，8,5,7,3であるから，長さ23（=8+5+7+3）のベクトルができることになる。

S1-2

25

S1-refa2.0~a2.3

S1-2

PE#0

PE#1

PE#2

PE#3

8
101.0
103.0
105.0
106.0
109.0
111.0
121.0
151.0

5
201.0
203.0
205.0
206.0
209.0

7
301.0
303.0
305.0
306.0
311.0
321.0
351.0

3
401.0
403.0
405.0

26

S1-ref

S1-2：局所⇒全体ベクトル生成：手順
PE#0

PE#1

stride(1)rcounts(1)displs(1)=0

stride(2)rcounts(2)PE#2

PE#(m-2) PE#(m-1)

stride(3)stride(m-1)

stride(m)rcounts(3)rcounts(m-1)

rcounts(m)displs(2)=displs(1) + stride(1)displs(m+1)=
displs(m) + stride(m)size(recvbuf)= displs(PETOT+1)= sum(stride)

S1-2

•
•
•
•

局所ベクトル情報を読み込む「rcounts」， 「displs」を作成する「recvbuf」を準備する
ALLGATHERV27

S1-ref

S1-2：局所⇒全体ベクトル生成（1/2）
s1-2.fimplicit REAL*8 (A-H,O-Z)
include 'mpif.h'
integer :: PETOT, my_rank, SOLVER_COMM, ierrreal(kind=8), dimension(:), allocatable :: VEC, VEC2, VECg
integer (kind=4), dimension(:), allocatable :: COUNT, COUNTindex
character(len=80):: filename
call MPI_INIT(ierr)
call MPI_COMM_SIZE (MPI_COMM_WORLD, PETOT, ierr )call MPI_COMM_RANK (MPI_COMM_WORLD, my_rank, ierr )if (my_rank.eq.0) filename= 'a2.0'
if (my_rank.eq.1) filename= 'a2.1'
if (my_rank.eq.2) filename= 'a2.2'
if (my_rank.eq.3) filename= 'a2.3'
open (21, file= filename, status= 'unknown')
read (21,*) N
allocate (VEC(N))
do i= 1, N
read (21,*) VEC(i)enddo

中身を書き出して見よう

allocate (COUNT(PETOT), COUNTindex(PETOT+1))
call MPI_allGATHER ( N, 1, MPI_INTEGER,&
COUNT, 1, MPI_INTEGER,&MPI_COMM_WORLD, ierr)
COUNTindex(1)= 0&
各PEにおけるベクトル長さの情報が&「COUNT」に入る（「rcounts」）

do ip= 1, PETOTCOUNTindex(ip+1)= COUNTindex(ip) + COUNT(ip)
enddo

S1-2

中身を書き出して見よう28

S1-ref

MPI_ALLGATHER

P#0 A0P#1 B0P#0 A0 B0 C0 D0
All gather

P#1 A0 B0 C0 D0P#2 C0P#2 A0 B0 C0 D0P#3 D0P#3 A0 B0 C0 D0• MPI_GATHER＋MPI_BCAST• call MPI_ALLGATHER (sendbuf, scount, sendtype, recvbuf,
rcount, recvtype, comm, ierr)
–
–
–
–
–
–
–
–

S1-2

sendbuf 任意
scount整数
sendtype 整数
recvbuf 任意
rcount整数
recvtype 整数
comm
整数
ierr
整数

I
I
I
O
I
I
I
O

送信バッファの先頭アドレス，
送信メッセージのサイズ
送信メッセージのデータタイプ
受信バッファの先頭アドレス，
受信メッセージのサイズ
受信メッセージのデータタイプ
コミュニケータを指定する
完了コード

29

S1-ref

S1-2：局所⇒全体ベクトル生成（2/2）
s1-2.f/c
do ip= 1, PETOTCOUNTindex(ip+1)= COUNTindex(ip) + COUNT(ip)
enddo「displs」に相当するものを生成。

allocate (VECg(COUNTindex(PETOT+1)))
VECg= 0.d0&&&

call MPI_allGATHERv( VEC , N, MPI_DOUBLE_PRECISION,
VECg, COUNT, COUNTindex, MPI_DOUBLE_PRECISION,MPI_COMM_WORLD, ierr)&&&

do i= 1, COUNTindex(PETOT+1)
write (*,'(2i8,f10.0)') my_rank, i, VECg(i)enddocall MPI_FINALIZE (ierr)
stop
end

PE#0

PE#1

stride(1)rcounts(1)displs(1)=0

S1-2

stride(2)rcounts(2)PE#2

PE#(m-2) PE#(m-1)

stride(3)stride(m-1)

stride(m)rcounts(3)rcounts(m-1)

rcounts(m)displs(2)=displs(1) + stride(1)size(recvbuf)= displs(PETOT+1)= sum(stride)

displs(m+1)=
displs(m) + stride(m)30S1-ref

S1-2：局所⇒全体ベクトル生成（2/2）
s1-2.f/c
do ip= 1, PETOTCOUNTindex(ip+1)= COUNTindex(ip) + COUNT(ip)enddoallocate (VECg(COUNTindex(PETOT+1)))
VECg= 0.d0&&&

「recvbuf」

call MPI_allGATHERv( VEC , N, MPI_DOUBLE_PRECISION,
VECg, COUNT, COUNTindex, MPI_DOUBLE_PRECISION,MPI_COMM_WORLD, ierr)&&&

do i= 1, COUNTindex(PETOT+1)
write (*,'(2i8,f10.0)') my_rank, i, VECg(i)enddocall MPI_FINALIZE (ierr)
stop
end

PE#0

PE#1

stride(1)rcounts(1)displs(1)=0

S1-2

stride(2)rcounts(2)PE#2

PE#(m-2) PE#(m-1)

stride(3)stride(m-1)

stride(m)rcounts(3)rcounts(m-1)

rcounts(m)displs(2)=displs(1) + stride(1)size(recvbuf)= displs(PETOT+1)= sum(stride)

displs(m+1)=
displs(m) + stride(m)31S1-ref

S1-2：局所⇒全体ベクトル生成（2/2）
s1-2.f/c
do ip= 1, PETOTCOUNTindex(ip+1)= COUNTindex(ip) + COUNT(ip)enddoallocate (VECg(COUNTindex(PETOT+1)))
VECg= 0.d0&&&

call MPI_allGATHERv( VEC , N, MPI_DOUBLE_PRECISION,
VECg, COUNT, COUNTindex, MPI_DOUBLE_PRECISION,MPI_COMM_WORLD, ierr)
do i= 1, COUNTindex(PETOT+1)
write (*,'(2i8,f10.0)') my_rank, i, VECg(i)enddocall MPI_FINALIZE (ierr)&&&call MPI_ALLGATHERV(sendbuf, scount, sendtype, recvbuf, rcounts, displs,
recvtype, comm, ierr)

stop
end

PE#0

PE#1

stride(1)rcounts(1)displs(1)=0

S1-2

stride(2)rcounts(2)PE#2

PE#(m-2) PE#(m-1)

stride(3)stride(m-1)

stride(m)rcounts(3)rcounts(m-1)

rcounts(m)displs(2)=displs(1) + stride(1)size(recvbuf)= displs(PETOT+1)= sum(stride)

displs(m+1)=
displs(m) + stride(m)32S1-ref

実 行（課題S1-2）
FORTRAN$ mpifrtpx –Kfast s1-2.f(modify “go4.sh”)$ pjsub go4.shC$ mpifccpx –Kfast s1-2.c(modify “go4.sh”)$ pjsub go4.sh

33

S1-ref

S1-2：結果
my_rank ID
0
1
0
2
0
3
0
4
0
5
0
6
0
7
0
8
0
9
0
10
0
11
0
12
0
13
0
14
0
15
0
16
0
17
0
18
0
19
0
20
0
21
0
22
0
23

S1-2

VAL
101.
103.
105.
106.
109.
111.
121.
151.
201.
203.
205.
206.
209.
301.
303.
305.
306.
311.
321.
351.
401.
403.
405.

my_rank ID
1
1
1
2
1
3
1
4
1
5
1
6
1
7
1
8
1
9
1
10
1
11
1
12
1
13
1
14
1
15
1
16
1
17
1
18
1
19
1
20
1
21
1
22
1
23

VAL
101.
103.
105.
106.
109.
111.
121.
151.
201.
203.
205.
206.
209.
301.
303.
305.
306.
311.
321.
351.
401.
403.
405.

my_rank ID
2
1
2
2
2
3
2
4
2
5
2
6
2
7
2
8
2
9
2
10
2
11
2
12
2
13
2
14
2
15
2
16
2
17
2
18
2
19
2
20
2
21
2
22
2
23

VAL
101.
103.
105.
106.
109.
111.
121.
151.
201.
203.
205.
206.
209.
301.
303.
305.
306.
311.
321.
351.
401.
403.
405.

my_rank ID
3
1
3
2
3
3
3
4
3
5
3
6
3
7
3
8
3
9
3
10
3
11
3
12
3
13
3
14
3
15
3
16
3
17
3
18
3
19
3
20
3
21
3
22
3
23

VAL
101.
103.
105.
106.
109.
111.
121.
151.
201.
203.
205.
206.
209.
301.
303.
305.
306.
311.
321.
351.
401.
403.
405.

34

S1-ref

S1-3：台形則による積分• 下記の数値積分の結果を台形公式によって求めるプログラ
ムを作成する。MPI_REDUCE，MPI_BCASTを使用して並列化を実施し，プロセッサ数を変化させた場合の計算時間を
測定する。

4
dx
2
0 1  x
1

S1-3

35

S1-ref

S1-3：台形則による積分プロセッサへの配分の手法
タイプA

タイプB

N
 を使うとすると必然的に
1 
2 f i x  f1  f N 1 2 i 2
 「タイプA」となるが・・・



S1-3

36

S1-ref

S1-3：台形則による計算
TYPE-A（1/2）：s1-3a.f
implicit REAL*8 (A-H,O-Z)
include 'mpif.h'
integer :: PETOT, my_rank, ierr, N
integer, dimension(:), allocatable :: INDEX
real (kind=8) :: dx
call MPI_INIT(ierr)
call MPI_COMM_SIZE (MPI_COMM_WORLD, PETOT, ierr )call MPI_COMM_RANK (MPI_COMM_WORLD, my_rank, ierr )allocate (INDEX(0:PETOT))INDEX= 0
if (my_rank.eq.0) thenopen (11, file='input.dat', status='unknown')
read (11,*) N
close (11)endif

“input.dat”で分割数Nを指定
中身を書き出して見よう:Ncall MPI_BCAST (N, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, ierr)
dx= 1.d0 / dfloat(N)nnn= N / PETOT
nr = N - PETOT * nnn
do ip= 1, PETOT
if (ip.le.nr) then
INEDX(ip)= nnn + 1
else
INDEX(ip)= nnn
endif
enddo

S1-3

各PEにおける小領域数が「nnn」（またはnnn+1）37

S1-ref

S1-3：台形則による計算
TYPE-A （2/2） ：s1-3a.f
do ip= 1, PETOTINDEX(ip)= INDEX(ip-1) + INDEX(ip)
enddo
Stime= MPI_WTIME()
SUM0= 0.d0
do i= INDEX(my_rank)+1, INDEX(my_rank+1)
X0= dfloat(i-1) * dx
X1= dfloat(i ) * dx
F0= 4.d0/(1.d0+X0*X0)
F1= 4.d0/(1.d0+X1*X1)
SUM0= SUM0 + 0.50d0 * ( F0 + F1 ) * dx
enddo
&

f0

call MPI_REDUCE (SUM0, SUM, 1, MPI_DOUBLE_PRECISION, MPI_SUM, 0,
MPI_COMM_WORLD, ierr)
Etime= MPI_WTIME()

f1

&

x0

x1

if (my_rank.eq.0) write (*,*) SUM, 4.d0*datan(1.d0), Etime-Stime
call MPI_FINALIZE (ierr)
stop
end

PE#0
INDEX(0)+1

S1-3

PE#1
INDEX(1)+1

PE#2
INDEX(2)+1

PE#(PETOT-1)
INDEX(3)+1

INDEX(PETOT-1)+1

INDEX(PETOT)=N
38

S1-ref

S1-3：台形則による計算
TYPE-B ：s1-3b.f
implicit REAL*8 (A-H,O-Z)
include 'mpif.h'
integer :: PETOT, my_rank, ierr, N
real (kind=8) :: dx
call MPI_INIT(ierr)
call MPI_COMM_SIZE (MPI_COMM_WORLD, PETOT, ierr )call MPI_COMM_RANK (MPI_COMM_WORLD, my_rank, ierr )if (my_rank.eq.0) thenopen (11, file='input.dat', status='unknown')
read (11,*) N
close (11)endifcall MPI_BCAST (N, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, ierr)
dx= 1.d0 / dfloat(N)Stime= MPI_WTIME()
SUM0= 0.d0
do i= my_rank+1, N, PETOT
X0= dfloat(i-1) * dx
X1= dfloat(i ) * dx
F0= 4.d0/(1.d0+X0*X0)
F1= 4.d0/(1.d0+X1*X1)
SUM0= SUM0 + 0.50d0 * ( F0 + F1 ) * dx
enddo&

call MPI_REDUCE (SUM0, SUM, 1, MPI_DOUBLE_PRECISION, MPI_SUM, 0,
MPI_COMM_WORLD, ierr)
Etime= MPI_WTIME()&

if (my_rank.eq.0) write (*,*) SUM, 4.d0*datan(1.d0), Etime-Stime

S1-3

call MPI_FINALIZE (ierr)
stop
end

39

S1-ref

コンパイル・実行（課題S1-3）
FORTRAN

タイプA$ mpifrtpx –Kfast s1-3a.f$ mpifrtpx –Kfast s1-3b.f(modify “go.sh”)$ pjsub go.sh

C

タイプB

$ mpifccpx –Kfast s1-3a.c$ mpifccpx –Kfast s1-3b.c(modify “go.sh”)$ pjsub go.sh40

S1-refgo.sh#!/bin/sh#PJM -L “node=1“#PJM -L “elapse=00:10:00“#PJM -L “rscgrp=school“#PJM -j#PJM -o “test.lst“
標準出力#PJM --mpi “proc=8“
mpiexec ./a.out

8分割“node=1““proc=8”16分割“node=1““proc=16”

32分割“node=2““proc=32”

64分割“node=4““proc=64”

192分割“node=12““proc=192”
41

S1-ref

S1-3：コンピュータにおける並列効果• ◆：N=106，●：108，▲：109，－：理想値• 1コアにおける計測結果（sec.）からそれぞれ算出• Strong Scaling– 全体問題規模固定– N倍のコアを使ってN分の1の計算時間で解けるのが理想• Weak Scaling– ノード当たり（コア当たり）問題規模固定– N倍のコアを使ってN倍の規模の問題
を同じ時間で解ける
S1-3
のが理想

42

S1-ref

理想値からのずれ• MPI通信そのものに要する時間– データを送付している時間– ノード間においては通信バンド幅によって決まる• Gigabit Ethernetでは 1Gbit/sec.（理想値）– 通信時間は送受信バッファのサイズに比例• MPIの立ち上がり時間– latency– 送受信バッファのサイズによらない• 呼び出し回数依存，プロセス数が増加すると増加する傾向– 通常，数～数十secのオーダー• MPIの同期のための時間– プロセス数が増加すると増加する傾向
43

S1-ref

理想値からのずれ（続き）• 計算時間が小さい場合（S1-３ではNが小さい場合）はこれらの効果を無視できない。
– 特に，送信メッセージ数が小さい場合は，「Latency」が効く。

44