Vol.-HPC- No. // 情報処理学会研究報告 integer :: array():[*] integer :: tmp() if (this_image() == ) then array(:)[] = tmp(:) tmp(:) = arrray(:)[] end if! Put co

京速コンピュータ「京」における

PGAS

モデルによる

気象コード

NICAM

の実装

中尾 昌広

1,2,a)

_{佐藤 三久}

1,2,3

概要：PGASモデル言語の1つであるXcalableMPを用いて全球雲解像モデルNICAMの通信モジュー

ルの実装を行い，その生産性と性能評価を京速コンピュータ「京」を用いて行った．生産性については，

MPIで記述されたNICAMの一対一通信をXcalableMPが提供するcoarray記法を用いて記述すること

で，その通信を簡易に表現することができた．また，「京」においてXcalableMPのcoarray記法で記述さ れたコードは，「京」が提供しているRDMA機能を利用するようにコード変換されるため，その通信の高 速化を見込むことができる．XcalableMPで実装を行ったNICAMを「京」を用いて性能評価した結果， 160ノード利用時に全体として約19%の高速化を達成することができた．

1.

はじめに

日本学術振興会・多国間国際研究協力事業[1]の研究プロ ジェクトの1つである「エクサスケール・コンピューティ ングによる精緻な気候シミュレーションの実現」では，各 国のスーパーコンピュータを用いた気候シミュレーション の研究が行われている．その研究の目的の1つに，エクサ スケール計算環境下におけるアプリケーション作成の生産 性向上がある． 分散メモリ環境におけるアプリケーションの作成には， 規模の大小を問わずMessage Passing Interface（MPI）[2]

が広く用いられている．しかしながら，近年，生産性の 面でMPIよりも有利なPartitioned Global Address Space

（PGAS）モデル[3]言語が普及しつつある．PGASモデル言 語の多くは片側通信と親和性が高く，またその片側通信は計 算ノードが持つRemote Direct Memory Access（RDMA）

機能を直接利用する場合がある．そのため，PGASモデル 言語を用いることで，MPIを用いた場合と比較して生産性 だけではなく性能も向上する可能性がある．

本稿では，PGASモデル言語の1つであるXcalableMP

（XMP）[4]を用いて，気候アプリケーションの1つである 全球雲解像モデルNICAM（Nonhydrostatic ICosahedral 1 _{筑波大学 計算科学研究センター}

Center for Computational Sciences, University of Tsukuba 2 _{理化学研究所 計算科学研究機構}

RIKEN Advanced Institute for Computational Science 3 _{筑波大学 大学院 システム情報工学研究科}

Graduate School of Systems and Information Engineering, University of Tsukuba a) _{[email protected]} Atmospheric Model）[5]の通信モジュールを実装し，その 生産性を評価する．また京速コンピュー タ「京」（以下， 京）を用いた性能評価も行う． 本稿の構成は次の通りである．2章ではPGASモデルの 概要について説明し，3章ではNICAMの概要について説 明する．4章ではXMPを用いたNICAMの実装について 述べ，5章では生産性と性能の評価を行う．6章で本稿の まとめと今後の課題について述べる．

2.

Partitioned Global Address Space モ

デル

PGASモデルは，各プロセスが透過的にアクセス可能な 領域を持つため，プロセス間通信を簡易に記述可能であ るという特徴がある．また，並列アプリケーションの作成 に必要なバリアやロックなどの機能も言語内に備えてい ることが多い．XMP以外のPGASモデル言語としては， SHMEM[6]，Global Arrays[7]，Coarray Fortran（CAF）[8]，

Titanium[9]，Unified Parallel C[10]，Chapel[11]，X10[12]

などがある．特にCAFはFortran 2008の標準規格に組み 込まれているため，今後の普及が期待できる． CAFでは，丸括弧を用いた通常の配列と，image番号 （CAFにおけるデータオブジェクトの識別番号）を指定す る角括弧の2つを用いてリモートのデータ領域をアクセス する．図1に，CAFのプログラミング例を示す．1行目で は，要素数10のinteger型のcoarrayを宣言している．5

行目では，image1がimage2のcoarrayに対してデータを 送信している（Put通信）．6行目では，image1がimage3

integer :: array(10):[*] integer :: tmp(10) if (this_image() == 1) then array(:)[2] = tmp(:) tmp(:) = arrray(:)[3] end if sync all 1 2 3 4 5 6 7 8 9 ! Define coarray ! Put communication ! Get communication ! Synchronization 図1 Coarray Fortranにおける通信の記述方法 Glevel-0 Glevel-1 図2 NICAMの格子点[5] のsync all文は，それより前に行われた片信通信を完了 させ，かつバリア同期を行う命令である．本稿では，図1 の5行目と6行目に示したような片側通信の記述方法を coarray記法と呼ぶ． XMPは，並列アプリケーションで広く用いられている FortranとC言語のそれぞれの言語拡張であり，XMPが 提供する指示文もしくはcoarray記法を用いて，通信を表 現する．また，Fortran版のXMPはCAFの上位互換とな るように設計されているため，CAFで記述されたプログ ラムをXMPとして動作させることが可能である．

3.

NICAM

NICAMは全球雲解像モデルの1つであり，Fortranと MPIライブラリで記述されている．NICAMが行う通信の 多くは一対一の隣接通信であるため，スケーラビリティが 非常に高いという特徴がある．本章では主に，並列計算に 必要な事柄について説明する． NICAMでは，全球に対して正二十面体格子を用いるこ とで，計算対象の点（水平格子点）を決定する．図 2に， NICAMの水平格子点の概念図を示す．まず，全球を正 二十面体に分割する．この状態をGlevel-0と呼ぶ（図 2 左）．その三角形のそれぞれの頂点が水平格子点である． Glevel-0のそれぞれの三角形を4分割した格子をGlevel-1 と呼ぶ（図2右）．このように再帰的に三角形を4分割し ていくことで，目的に応じた解像度をユーザが設定するこ とができる．再帰回数がnの場合の格子をGlevel-nと呼 ぶ．Glevel-nの場合の水平格子点数は10∗ 4n_{+ 2の式で} 計算することができる． 次に，並列化を行う場合の各プロセスが担当する領域の 設定方法について説明する．まず，Glevelと同様に，全球

Glevel-3, Rlevel-0 Glevel-3, Rlevel-1

図3 NICAMの領域分割（領域の一部のみ，太枠で囲んでいる）[5] を正二十面体に分割する．次に，その正二十面体の三角形 を2つずつ合わせて四角形の領域を作成する．この領域 が，各プロセスが担当する領域である．そして，Glevelと 同様に，再帰的に領域を4分割していくことで，多くの領 域を設定することが可能である．最初の分割された領域を Rlevel-0と呼び，次に4分割された領域をRlevel-1と呼ぶ． Rlevel-0の場合の領域数は10であり，Rlevel-1の場合の領 域数は40である．図3に例を示す．GlevelとRlevelは個 別に設定することが可能である．Rlevel-nの場合の領域数 は10∗ 4n_{の式で計算することができる．並列計算に用い} るプロセス数は，Rlevelによって設定された領域数の約数 である必要がある．例えば，Rlevel-0の場合の領域数は10 であるため，ユーザが利用できるプロセス数は1，2，5，10 となる．

4.

XcalableMP による NICAM の実装

4.1 関連研究 NICAMは地球シミュレータを用いて開発が行われてき たため，NICAMのコードはベクトル計算機用に最適化さ れている．そのため，京の性能を引き出すためのNICAM に対する最適化作業が現在進められている．[13]では，地 球シミュレータで動作していたコードに対してキャッシュ 最適化などを行うことにより，京上におけるNICAMの性 能効率を2倍以上に高めている．また[14]では，京の3次 元トーラスネットワークに対して，通信のホップ数が少な くなるようなプロセスの割り当て手法の提案が行われて いる． 4.2 最適化の方針 NICAMの一対一通信をXMPのcoarray記法による片 側通信によって記述することで，コードの簡易化を図る． また，XMPの実装において，京のRDMA機能を用いて片 側通信を実行させることにより，高速化も図る．表 1に， 京が提供している拡張RDMAインタフェースの一覧を示 す[15]．このRDMAインタフェースはC言語の関数とし て定義されている．C言語で実装されたXMPのランタイ ムライブラリから，表1の各関数を呼び出すことによって， 京のRDMA機能を直接用いることができる．

表1 拡張RDMAインタフェース[15]

関数名 機能

FJMPI Rdma init 拡張RDMAインターフェースの初期化

FJMPI Rdma finalize 拡張RDMAインターフェースの終了処理

FJMPI Rdma reg mem メモリの登録

FJMPI Rdma dereg mem メモリの登録解除

FJMPI Rdma get remote addr リモートDMAアドレスの獲得

FJMPI Rdma put RDMA WRITE通信

FJMPI Rdma get RDMA READ通信

FJMPI Rdma poll cq RDMA完了確認

MPI Isend/Irecv RDMA 10 10 10 10 10 10 10 10 5.0 4.5 4.0 3.5 3.0 2.5 2.0 1.5 1.0 0.5 0 Ba n d w id th (G Byt e /s)

Transfer Size (Byte)

0 1 2 3 4 5 6 7 図4 RDMAとMPI関数との通信性能比較 10 10 10 10 10 10 10 10 3.0 2.5 2.0 1.5 1.0 0.5 0 Pe rf o rma n ce R a ti o

Transfer Size (Byte)

0 1 2 3 4 5 6 7

RDMA Bandwidth/ MPI_Isend/Irecv Bandwidth

図5 RDMAとMPI関数との通信性能比

4.3 予備実験

NICAMの一対一通信は，MPI Isend/Irecv関数を用い

て実装されている．そこで，まず京のRDMA機能を用い た通信とMPI Isend/Irecv関数との性能比較を行う． 2ノード間でpingpongを行うプログラムを京のRDMA 機能を用いた通信およびMPI Isend/Irecv関数を用いて それぞれ作成し，京において実行した．結果を図 4に示 す．また，図4における性能比を示すグラフを図5に示す． 図5の縦軸の値が1.0以上であれば，RDMA機能を用いた 通信の方の性能が高い．図4と図5より，すべての転送サ イズにおいてRDMA機能を用いた方がMPI Isend/Irecv

関数を用いた場合よりも高速であることがわかる．また

do i=1, recv_num

call mpi_irecv(recvbuf(1,i), recv_count(i),

mpi_double_precision, sourcerank(i), ...) enddo

...

do i=1, send_num

call mpi_isend(sendbuf(1,i), send_count(i), mpi_double_precision, destrank(i) ...) enddo ... call mpi_waitall(...) 1 2 3 4 5 6 7 8 9 10 11 MPI (Original)

real(8) :: recvbuf(maxdatasize_r, romax(halomax)):[*] real(8) :: sendbuf(maxdatasize_s, somax(halomax)):[*]

...

! Obtain information of destination position -> dest_position() ... do i=1, send_num recvbuf(1:send_count(i), dest_position(i))[destrank(i)] = sendbuf(1:send_count(i), i) enddo ... sync all 1 2 3 4 5 6 7 8 9 10 11 XcalableMP 図6 XcalableMPによるNICAMの通信モジュールの実装例（変 数名などは一部変更） 図5より，転送サイズが105_{Byte以下の場合，RDMA機} 能を用いた通信は特に高速であることがわかる． 4.4 XcalableMPによるNICAMの通信モジュールの 実装

NICAMはFortranで記述されているため，Fortran版の

XMPを用いてNICAMの通信モジュールの実装を行う． 図6に，オリジナルのNICAMからFortran版のXMPに 変換したコードの一部を示す（図における分かりやすさを 優先させたため，実際のコードとは変数名や配列の次元数 は異なる）． 図6下のXMPのコードにおいて，1行目と2行目では 配列recvbufとsendbufをCoarrayとして宣言している．

XMP（およびCAF）の文法上は，配列sendbufはCoarray

である必要はないが，京では送信用バッファのアドレスも

RDMAを用いる領域として登録する必要があるため，配列 sendbufもCoarrayとして宣言している．4行目（実際は複

数行．紙面の都合上省略している）では，送信するデータ の格納先の情報を，MPI Allgather通信を用いることで各 プロセスから取得し，配列dest positionに保存している． この操作は，オリジナルのMPI版では必要のない操作で ある．なお，NICAMでは通信を行う対象プロセスとその データの格納先のアドレスはプログラム実行中は不変であ るため，この操作はプログラム開始時に一度行うだけで良 い．6行目と7行目では，Put通信を行っている．最後に 11行目では，片側通信の完了処理とバリア同期を行って いる． XMPで実装したNICAMとオリジナルのNICAMとの 他の違いには，配列recvbufとsendbufの要素数が異なる 点が挙げられる．オリジナルのNICAMでは，プロセス毎 に配列recvbufとsendbufの要素数は異なる．XMP（およ びCAF）の規格上，各プロセスで宣言するCoarrayは同 じ要素数である必要がある．そのため，今回の実装では， 各プロセスで最大となるそれぞれの配列の要素数を用い てCoarrayの宣言を行っている．ただし，配列recvbufと sendbufの要素数は各プロセスでほぼ同じであるため，最大 要素数を用いることによる利用メモリ量の増加はデメリッ トにはならない．

5.

評価

5.1 生産性 生産性を表す指標には行数がよく用いられる．今回の実 装において行数に着目した場合，通信のコードが全体を占 める割合は僅かであるため，全体としての行数はほとんど 変化しない．しかしながら，coarray記法によるXMPの実 装では，通信の記述は送信元プロセスのみでよく，またオ リジナルのNICAMで用いられているMPI関数の呼び出 しよりも直感的に通信の記述が行える．以上のことから， XMPのcoarray記法を用いた実装の方が，より簡潔に通 信を表現できていると言える． 5.2 性能 京を用いてXMPで実装したNICAMとオリジナルの NICAMとの性能比較を行う．ただし，Fortran版のXMP の片側通信機能は未実装であるため，今回の性能評価にお いては，図6下に示したXMPのコードからXMPコンパ イラが変換すると考えられる表1のRDMA関数を直接呼 び出す形で実装を行っている．京のシステム概要を表2に 示す． 実験で用いたNICAMの並列計算に関係するパラメータ として，Glevelは5（格子点数は10242）に固定し，Rlevel は0，1，2と変化させて実験を行った．並列数は各Rlevel の最大並列数である10，40，160に設定した．NICAMは OpenMPおよびコンパイラの自動並列化機能によりスレッ ド並列化されて実行されるため，京の各計算ノードに対し Rlevel-0 (10 nodes) Rlevel-1 (40 nodes) Rlevel-2 (160 nodes) 8.0 7.0 6.0 5.0 4.0 3.0 2.0 1.0 0 T ime (s) (0.46)3.60 (0.56)3.72 7.53 (0.41) 7.60 (0.31) 2.35 (0.93) RDMAăઌ̲࢘в RDMAăઌ࢘ MPIăઌ̲࢘в MPIăઌ࢘ 2.79 (1.22) 図7 性能評価と内訳（括弧内の数値は通信に要する時間） て1プロセスを割り当て，各プロセスは8スレッドで動作 するように設定した． 他のパラメータとして，シミュレーションを行う時間変 化量とステップ数がある．今回の実験では時間変化量は 1200秒，ステップ数は12に設定した．すなわち，14400秒 （4時間）のシミュレーションである．またNICAMではス テップ毎に中間ファイルを生成し，さらに最終ステップに リスタート用のファイルを生成する．今回の実験では，低 いGlevelを用いたため，これらのファイルの生成に要する 時間が全体の計算時間に与える影響が大きい．そのため， これらのファイルの生成は行わないようにした． 性能評価の結果を図 7に示す．図7の括弧内の数値は， 今回実装を行った通信に関係する箇所の時間（通信の開始 から同期をとるまでの時間．その間に計算も行われている） である．図 7より，利用するノード数が多いほど，XMP で実装したNICAMの方がオリジナルのNICAMよりも性 能がより高くなることがわかる．160ノードを利用した際 のXMPの実装とオリジナルとの速度差は，通信部分だけ では約31%の速度向上，全体としては約19%の速度向上で あった． NICAMでデータ交換を行う箇所は各領域の袖の部分で あり，Rlevelが1つ上がる毎に，各領域に存在する水平格 子数は1/4になる一方，データ通信量は1/2にしかならな い．すなわち，Rlevelが上がる毎に通信時間の占める割合 が増えるため，Rlevelが大きい方がより通信の性能差が現 れやすくなる． Rlevel-0 に お け る1回 の 通 信 の 平 均 転 送 サ イ ズ は 約

23KByteであり，Rlevel-1の平均転送サイズは約9KByte，

Rlevel-2の平均転送サイズは約5KByteである．図5より， それぞれの転送サイズにおけるRDMAとMPIとの速度差 は約2.5倍である．図7の括弧内のそれぞれの数値が2.5 倍も差はない理由は，その数値は同期待ちの時間が含まれ ており，また通信と同時に計算も行われているからと考え られる．

表2 京のシステム概要

CPU SPARC64 VIIIfx 2.0GHz, 8Cores/Socket, 128GFlops

Memory DDR3 SDRAM 16GB, 64GB/s/Socket

Network Torus fusion six-dimensional mesh/torus network, 5GB/s Compiler Fujitsu Fortran Compiler Version K-1.2.0-13

Communication Library Fujitsu MPI Version K-1.2.0-13

6.

まとめと今後の課題

本稿では，PGASモデル言語の1つであるXMPを用い て，NICAMの通信モジュールの実装を行い，京の上で性 能評価を行った．MPI関数で記述されている一対一通信 を，XMPのcoarray記法を用いることによって，京の持 つRDMA機能を直接利用でき，通信の高速化を図ること ができる．京の計算ノードを最大160ノード用いて性能評 価を行った結果，全体として約19%の高速化を達成するこ とができた．また，一対一通信にcoarray記法を用いるこ とで，ソースコードの簡易化を行うことができた． 今後の課題として，Fortran版のXMPのCoarray実装 を進めることと，より解像度が高いデータとより多くの計 算ノードを用いて性能評価を行うことが挙げられる．また， 図6における今回の実装では全プロセス間でバリア同期を 行うsync all文を用いたが，実際は通信を行う相手間の 同期のみで計算を進めることができる．すなわち，部分的 なバリア同期を行うことにより，性能向上を図ることがで きる．部分的なバリア同期を行うための命令として，XMP

やCAFではsync images文が規格されている．さらに， [14]のような京のネットワークに最適なプロセスマッピン グ方法と本実装とを併用することで，さらなる高速化が可 能になると考えられる． 謝辞 本研究を遂行するにあたり，NICAMを提供して くださった東京大学の佐藤正樹先生を始めとするNICAM 開発グループの皆様に感謝の意を表します．また，NICAM の利用方法，チューニングなどについてアドバイスを下 さった筑波大学計算科学研究センターの寺崎康児研究員に 感謝の意を表します． 本研究は，日本学術振興会・多国間国際研究協力事業 「エクサスケール・コンピューティングによる精緻な気候 シミュレーションの実現」の支援によって行われました． 参考文献 [1] 日 本 学 術 振 興 会 多 国 間 国 際 研 究 協 力 事 業:事 業 概 要 ． http://www.jsps.go.jp/j-bottom/01 b gaiyo.html [2] M. Snir, S. Otto, S. Huss-Lederman, D. Walker, and J.

Dongarra, MPI-The Complete Reference, Volume 1: The MPI Core, 2nd ed. Cambridge, MA, USA: MIT Press, 1998.

[3] PGAS - Partitioned Global Address Space Languages．

http://www.pgas.org [4] http://www.xcalablemp.org

[5] Satoh, M., T. Matsuno, H. Tomita, H. Miura, T. Nasuno,

S. Iga (2008), “Nonhydrostatic Icosahedral Atmospheric Model (NICAM) for global cloud resolving simulations.” Journal of Computational Physics, the special issue on Predicting Weather, Climate and Extreme events, 227, 3486-3514, doi:10.1016/j.jcp.2007.02.006.

[6] B. Chapman, T. Curtis, S. Pophale, S. Poole, J. Kuehn, C. Koelbel, and L. Smith, “Introducing openshmem: Shmem for the pgas community,” in Proceedings of the Fourth Conference on Partitioned Global Address Space Programming Model, ser. PGAS ’10. New York, NY, USA: ACM, 2010, pp. 2:1 – 2:3.

[7] J. Nieplocha, R. J. Harrison, and R. J. Littlefield, “Global arrays: A non-uniform-memory-access program-ming model for high-performance computers,” THE JOURNAL OF SUPERCOMPUTING, vol. 10, pp. 169 – 189, 1996.

[8] R. W. Numrich and J. Reid, “Co-array fortran for par-allel programming,” SIGPLAN Fortran Forum, vol. 17, no. 2, pp. 1 – 31, Aug. 1998.

[9] K. Yelick, L. Semenzato, G. Pike, C. Miyamoto, B. Li-blit, A. Krish- namurthy, P. Hilfinger, S. Graham, D. Gay, P. Colella, and A. Aiken, “Titanium: A high-performance Java dialect,” in ACM 1998 Workshop on Java for High-Performance Network Computing. New York, NY 10036, USA: ACM Press, 1998.

[10] U. Consortium, “UPC Language Specifications,” Berke-ley National Laboratory, Tech. Rep. LBNL-59208, 2005. [11] B. Chamberlain, D. Callahan, and H. Zima, “Parallel programmability and the chapel language,” Int. J. High Perform. Comput. Appl., vol. 21, no. 3, pp. 291–312, Aug. 2007.

[12] V. Saraswat, B. Bloom, I. Peshansky, O. Tardieu, and D. Grove, “X10 language specification,” 2013, http://x10.sourceforge.net/documentation/ languagespec/x10-231.pdf. [13] 八代尚，“全球非静力学モデルNICAMの現在と今後の開 発計画について”，地球流体データ解析・数値計算ワーク ショップ，http://www.gfd-dennou.org/arch/davis/ workshop/2012-12-12/yashiro 20121212.pdf， 2012年 [14] 小玉知央，寺井優晃，野田暁，山田洋平，佐藤正樹，清木 達也，伊賀晋一，富田浩文，南一生，“正二十面体格子にお けるノードマッピング手法の開発と評価”，京コンピュー タ・シンポジウム2012，2012年

[15] Parallelnavi for MP10 V1.0，Parallelnavi Technical Com-puting Language MPI使用手引書，2013年