XMPによる並列化実装 2
中尾 昌広（理化学研究所 計算科学研究機構）

KOBE HPC Spring School 2018

目次
グローバルビュープログラミング
ローカルビュープログラミング（時間があれば）2用いるコード./xmp-sample/global-view

Exercise 1

Exercise 2C

FortranSerialinit.cinit.f90

XMP

xmp_init.c

xmp_init.f90

Seriallaplace.c

laplace.f90

XMP

xmp_laplace.c

xmp_laplace.f90

3

課題1
分散配列a[]を定義し，ループ文を2ノードで並列実行させよ[C] init.c[F] init.f90#include <stdio.h>

program init

int a[10];integer :: a(10), i

int main(){

do i=1,10a(i)=ifor(int i=0;i<10;i++)end doa[i] = i+1;
do i=1,10for(int i=0;i<10;i++)printf("%d¥n", a[i]);print *, a(i)end doreturn 0;

end program init

}

4

課題1
下記のファイルは部分的に並列化しています
XMP指示文をxmp̲init.c or xmp̲init.f90に挿入してください[C] xmp̲init.c[F] xmp̲init.f90#pragma xmp nodes p[2]#pragma xmp template t[10]program init!$xmp nodes p(2)#pragma xmp distribute t[block] onto p
int a[10];!$xmp template t(10)!$xmp distribute t(block) onto p[align directive]

integer :: a(10), i

int main(){[loop directive][align directive]for(int i=0;i<10;i++)a[i] = i+1;[loop directive]
do i=1,10a(i)=i[loop directive]for(int i=0;i<10;i++)end do[loop directive]

printf(”[%d] %d¥n", xmpc_node_num(), a[i]);return 0;

do i=1,10
print *, xmp_node_num(), a(i)}end doend program init

5

実行手順
コンパイル
[C] $ xmpcc xmp̲init̲c -o xmp̲init[F] $ xmpf90 xmp̲init.f90 -o xmp̲init
2ノードで実行6課題2：ラプラス関数
初期値は0
初期値は下記の通り

赤のセルは上下左右のセルの値を使って更新される．
7

逐次コードまずはlaplace.c or laplace.f90をエディタで開いてください
下記のようにコンパイルしてください[C] $ gcc laplace.c -lm[F] $ gfortran laplace.f90
フロントエンドノードで実行$ ./a.out
Veriﬁcation = 5.54885...

8並列化のパターン[Pattern 1] 両辺で配列のインデックスが同じ場合[C] u[j][i] = uu[j][i];[F] u(i,j) = u(i,j)
[C]

[F]

#pragma xmp loop (j,i) on t[j][i]!$xmp loop (i,j) on t(i,j)

for(int j=0;j<10;j++)for(int i=0;i<10;i++)do j=1,10
do i=1,10u[j][i] = uu[j][i];u(i,j) = u(i,j)
end doend do9並列化のパターン[Pattern 2] 隣接のセルを参照する場合[C] u[j][i] = uu[j-1][i] + uu[j][i-1] + ...;
[F] u(i,j) = uu(i,j-1) + uu(i-1,j) + ...
shadow指示文を使って，袖領域を分散配列に追加する
reﬂect指示文を使って，ループ文の前で袖領域の更新を行う
[C]

[F]

#pragma xmp shadow uu[1:1][1:1]:!$xmp shadow uu(1:1,1:1):#pragma xmp reflect (uu)!$xmp reflect (uu):#pragma xmp loop (j,i) on t[j][i]:!$xmp loop (i,j) on t(i,j)

for(int j=0;j<10;j++)for(int i=0;i<10;i++)do j=1,10
do i=1,10u[j][i] = uu[j-1][i] + uu[j][i-1] + ...

u(i,j) = u(i,j-1) + uu(i-1,j) +end doend do10並列化のパターン[Pattern 3] reductionが必要な場合[C] s += abs(uu[j][i]-u[j][i]);[F] s = s + abs(uu(i,j)-u(i.j))[C]#pragma xmp loop (j,i) on t[j][i] reduction (+:s)
for(int j=0;j<10;j++)for(int i=0;i<10;i++)s += abs(uu[j][i] - u[j][i]);[F]!$xmp loop (i,j) on t(i,j) reduction(+:s)
do j=1,10
do i=1,10s = s + abs(uu(i,j) - u(i,j))
end doend do11

リモートデータの参照
逐次実行

u[j][i] = (u[j+1][i]+u[j-1][i]+u[j][i+1]+u[j][i-1])/4;

並列実行

参照に通信が必要

12

リモートデータの参照
1. shadow指示文は袖領域を分散配列に追加する
2. reﬂect指示文は，隣接ノードから分散配列の端の値を，自ノードの袖領域にコピーする
3. loop指示文は，隣接セルの参照を含む
for文を処理する13

課題2
xmp̲laplace.c or xmp̲laplace.f90を元に，XMPで並列化せよ
2次元ノード集合と2次元テンプレートを用いる
4ノードと8ノードで実行させる
逐次バージョンの結果と比べる

14

XMP version Laplace equation[loop directive]
pattern 1for(j = 1; j < N2-1; j++)for(i = 1; i < N1-1; i++)u[j][i] = sin((double)i/N1*M_PI)+ cos((double)j/N2*M_PI);#include <stdio.h>#include <stdlib.h>#include <math.h>#define N1 64#define N2 64
double u[N2][N1],uu[N2][N1];for(k = 0; k < niter; k++){/* old <- new */[loop directive]for(j = 1; j < N2-1; j++)for(i = 1; i < N1-1; i++)uu[j][i] = u[j][i];#pragma xmp nodes p[*][4]#pragma xmp template t[N2][N1][distribute directive][align directive][align directive][shadow directive]
int main(int argc, char **argv){
int j,i,k,niter = 100;
double value = 0.0;
#pragma xmp loop (j,i) on t[j][i]for(j = 0; j < N2; j++){for(i = 0; i < N1; i++){u[j][i] = 0.0;
uu[j][i] = 0.0;}}

pattern 1

pattern 1[reflect directive]
pattern 2[loop directive]for(j = 1; j < N2-1; j++)for(i = 1; i < N1-1; i++)u[j][i] = (uu[j-1][i] + uu[j+1][i] +uu[j][i-1] + uu[j][i+1])/4.0;}/* check value */value = 0.0;
#pragma xmp loop (j,i) on t[j][i] [reduction clause]for(j = 1; j < N2-1; j++)for(i = 1; i < N1-1; i++)value += fabs(uu[j][i]-u[j][i]);pattern 3

15

目次
グローバルビュープログラミング
ローカルビュープログラミング（時間があれば）

16

ローカルビュープログラミング
Coarray記法について
課題1 : スカラー値のPut/Get（実行のみ）
課題2 : 配列のPut/Get（実行のみ）
課題3 : 行列積の作成

17

Coarray記法
片側通信 (Put/Get)array[start:length[:step]]
ノード1が持っている b[3], b[4], b[5]の値が，ノード番号0のa[0], a[1], a[2]にコピーされる

image 0image 1if(xmpc_this_image() == 0)a[0:3] = b[3:3]:[1]// Get

a[10]

a[10]

角括弧はノード番号を示す

b[10]

b[10]

int a[10];int b[10]:[*]; // b is a coarray[C]:（Fortran is 1-origin, C is 0-origin）

3

5
18

Array section in XMP/C (1/2)
シンタックス
array[base : length :step]base, length, stepはint型
Base:省略された場合は0になる
Length :省略された場合は，残りの要素数になる
Step :
省略された場合は1になる

19

Array section in XMP/C (2/2)
int A[100]:[*]のCoarrayが宣言されている場合

A[10:10]

A[10] から A[19]の10要素

A[10:]

A[10] から A[99]の90要素

A[:10]

A[0] から A[9]の10要素

A[10:5:2]

A[10], A[12], A[14], A[16], A[18]の5要素

A[:]

Aの全要素 (A[0] から A[99])20

Coarray記法
片側通信 (Put/Get)

ノード2が持っているb(3:5)はノード1のa(1:3)にコピーされる

image 1

image 2if(this_image() == 1) thena(1:3) = b(3:5)[2]// Get

a(10)

a(10)

角括弧はノード番号を示す

b(10)

b[10]

integer :: a(10)integer :: b(10)[*] // b is a coarray[F]:（Fortran is 1-origin, C is 0-origin）

3

5
21

同期
[C]

void xmp_sync_all(int *status)[F]

sync allimage 0image 1

Put
sync_allsync_all

これまでに実行した片側通信が完了させ，かつノード間でバリア同期を行う（MPI̲Win̲fence+ MPI̲Barrierに相当）

22

課題1
スカラ変数に対するPut/Get[C] xmpcc coarray̲scalar.c -o coarray̲scalar.x[F] xmpf90 coarray̲scalar.f90 -o coarray̲scalar.x
2プロセスで実行せよ
[C]

[F]

if(xmpc_this_image() == 0){
tmp = val:[1]; // Get
val:[1] = val; // Put}
xmp_sync_all(NULL);

if(this_image() == 1) thentmp = val[2] ! Getval[2] = val ! Put
end ifsync all23

課題1 : 結果 (XMP/C)
xmp_sync_all(NULL);
if(xmpc_this_image() == 0){
tmp = val:[1]; // Get
val:[1] = val; // Put}
xmp_sync_all(NULL);

image 0val = 1
tmp = 0val = 2
tmp = 0xmp_sync_all

Get val:[1] to tmp

結果tmp = 2

[START] My image is 0, val = 1 tmp = 0[START] My image is 1, val = 2 tmp = 0

Put valto val:[1]

[END]
[END]

image 1

My image is 0, val = 1 tmp = 2
My image is 1, val = 1 tmp = 0

val = 1
xmp_sync_all24課題1：結果 (XMP/Fortran)
sync allif(this_image() == 1) thentmp = val:[2]; // Get
val:[2] = val; // Put
end ifsync allimage 1
val = 1
tmp = 0val = 2
tmp = 0sync allGet val:[2] to tmp

tmp = 2

結果
[START] My image is 1, val = 1 tmp = 0[START] My image is 2, val = 2 tmp = 0
[END]
[END]

image 2

My image is 1, val = 1 tmp = 2
My image is 2, val = 1 tmp = 0

Put val toval:[2]

val = 1sync all25

課題2
配列に対するPut/Get[C] xmpcc coarray̲vector.c -o coarray̲vector.x[F] xmpf90 coarray̲vector.f90 -o coarray̲vector.x
2プロセスで実行せよ[C] ノード 0, [F] ノード 1[C] ノード 1, [F] ノード 2

a[10] = {0, 1, .., 9};b[10] = {0, 1, .., 9};c[10][10] = {{0, 1, .., 9},{10, 11, .., 19},...
{90, 91, .., 99}};a[10] = {10, 11, .., 19};b[10] = {10, 11, .., 19};c[10][10] = {{100, 101, .., 109},{110, 111, .., 119},...
{190, 191, .., 199}};26

課題 2 : 結果
[C]

[C]

if(xmpc_this_image() == 0){a[0:3] = a[5:3]:[1]; // Get}

a[0] = 15a[1] = 16a[2] = 17a[3] = 3a[4] = 4a[5] = 5a[6] = 6a[7] = 7a[8] = 8a[9] = 9[F]if(this_image() == 1) thena(1:3) = a(6:8)[2] ! Get
end if[F]a(1) = 15a(2) = 16a(3) = 17
a(4) = 3a(5) = 4a(6) = 5a(7) = 6a(8) = 7a(9) = 8a(10) = 927

課題 2 : 結果
[C]

[C]

[F]

if(xmpc_this_image() == 0){b[0:5:2] = b[0:5:2]:[1];// Get}b[0] = 10b[1] = 1b[2] = 12b[3] = 3b[4] = 14b[5] = 5b[6] = 16b(1) = 10b(2) = 1b(3) = 12b(4) = 3b(5) = 14[F]if(this_image() == 1) thenb(1:10:2) = b(1:10:2)[2]
end if! Getb[7] = 7b[8] = 18
b[9] = 9b(6) = 5b(7) = 16b(8) = 7b(9) = 18
b(10) = 928

課題 2 : 結果
[C]

[F]

if(xmpc_this_image() == 0){c[0:5][0:5]:[1] = c[0:5][0:5]; // Put}if(this_image() == 1) thenc(1:5,1:5)[2] = c(1:5,1:5) // Put
end if0 1 2 3 4 105 106 107 108 109
10 11 12 13 14 115 116 117 118 119
20 21 22 23 24 125 126 127 128 129
30 31 32 33 34 135 136 137 138 139
40 41 42 43 44 145 146 147 148 149
150 151 152 153 154 155 156 157 158 159
160 161 162 163 164 165 166 167 168 169
170 171 172 173 174 175 176 177 178 179
180 181 182 183 184 185 186 187 188 189
190 191 192 193 194 195 196 197 198 199
29

課題 3
行列積の並列化
C=AxB

並列化の手順4ノードで実行する

for(i=0;i<N;i++)for(j=0;j<N;j++)for(k=0;k<N;k++)c[i][k] += a[i][j] * b[j][k];まず，ノード 0 in C（ノード1 in Fortran）はa[][]とb[][]のすべてのデータを代入する
次に，ノード 0 in C（ノード1 in Fortran）はa[][]とb[][]の一部を他のノードにPutする
そして，4ノードは a[][]とb[][]を使ってc[][]の計算を行う
最後に，ノード 0 in C（ノード1 in Fortran）は，他のノードからc[][]の一部を自分のc[][]に集める

30

課題 3
Calculation by sub-matrix

C00

C01

A00

A01

=
C10

C11

B00

B01

B10

B11

×
A10

A11

XMP/C

XMP/Fortran

C00 = A00 × B00 + A01 × B10

ノード 0, ノード 1

C01 = A00 × B01 + A01 × B11
C10 = A10 × B00 + A11 × B10

ノード 1, ノード 2

C11 = A10 × B01 + A11 × B11

ノード 3, ノード 4

ノード 2, ノード 3

31

課題 3
Calculation by sub-matrix

C00

C01

A00

A01

=
C10

C11

B00

B01

B10

B11

×
A10

A11

XMP/C

XMP/Fortran

C00 = A00 × B00 + A01 × B10

ノード 0, ノード 1

C01 = A00 × B01 + A01 × B11
C10 = A10 × B00 + A11 × B10

ノード 1, ノード 2

C11 = A10 × B01 + A11 × B11

ノード 3, ノード 4

ノード 2, ノード 3

32

課題 3
matmul.c or matmul.f90を元に，XMPで並列化を行う．
下記の赤字の箇所だけ，実装を行う．
1. (To compare time and veriﬁcation) serial matrix multiply
2. In the function init̲dmat(), Initialize arrays
3. In the function move̲data(),
The root image puts sub-matrix A[][] and B[][] to another images[C] Image 1 requires A[0:N/2][0:N] and B[0:N][N/2:N/2][F] Image 2 requires A(1:N,1:N/2) and B(N/2+1:N,1:N)
4. In the function mul̲dmat(), execute a matrix multiply (C = A x B)
5. In the function gather̲data(), gather results from each node
6. In the function verify(), perform veriﬁcation
If the value is the same result as the sequential version, the result is OK
Please compare the measurement time33