応用数理工学特論 線形計算と ハイパフォーマンスコンピューティング 応用数理工学特論 線形計算と ハイパフォーマンスコンピューティング 第3回 計算理工学専攻　張研究室 山本有作 日立製作所の山本有作です。 「～」について発表いたします。

前回の概要 「並列計算機による高性能計算」 １. 並列計算機の種類と特徴 ２. 共有メモリ型並列計算機 プログラミングモデル 高性能化の技法 １.　並列計算機の種類と特徴 ２.　共有メモリ型並列計算機 プログラミングモデル 高性能化の技法 本発表では，はじめに，研究の背景を述べてから，スパースソルバの概要，並列化手法，そして本研究で工夫した点の一つであるRISCプロセッサ向けの最適化についてご説明します。最後に，並列計算機SR2201上での性能評価とまとめを述べます。

今回の概要 「並列計算機による高性能計算」 ３. 分散メモリ型並列計算機 「連立一次方程式の高性能解法 （密行列の場合）」 １. LU分解 ３.　分散メモリ型並列計算機 プログラミングモデル 高性能化の技法 「連立一次方程式の高性能解法　（密行列の場合）」 １.　LU分解 ２.　LU分解のブロック化 ３.　LU分解の並列化 本発表では，はじめに，研究の背景を述べてから，スパースソルバの概要，並列化手法，そして本研究で工夫した点の一つであるRISCプロセッサ向けの最適化についてご説明します。最後に，並列計算機SR2201上での性能評価とまとめを述べます。

3.1　分散メモリ型並列計算機のプログラミング 通信ライブラリMPI MPIのプログラミングモデル MPIの関数 行列乗算の例

通信ライブラリMPI MPI（Message Passing Interface） MPIの機能 プロセッサ間通信を行うためのサブルーチンの集合 MPIの機能 １対１通信 ブロードキャスト 総和演算・最大値 その他

MPIのプログラミングモデル SPMD（Single Program Multiple Data） 分散メモリ 全プロセッサが共通のプログラムを実行 処理するデータはプロセッサ毎に異なる。 各プロセッサは固有の番号（0, 1, … , p-1）を持ち，プログラム中で自分の番号を参照して処理を行う。 分散メモリ 各プロセッサは自分の持つローカルメモリのみをアクセス可能 他プロセッサのメモリ上にあるデータが必要な場合は，プロセッサ間通信により取得

MPIの関数 起動／終了 環境情報の取得 １対１の送受信 集団通信 MPI_INIT ： MPIの初期化 MPI_FINALIZE ：　MPIの終了 環境情報の取得 MPI_COMM_SIZE ：　全プロセッサ台数の取得 MPI_COMM_RANK ：　プロセッサ番号の取得 １対１の送受信 MPI_SEND ：　データの送信 MPI_RECV ：　データの受信 集団通信 MPI_BCAST ：　データのブロードキャスト MPI_REDUCE ：　リダクション演算（総和／最大値など）

行列乗算の例 (1) × = 問題 PU0 PU3 PU0 PU1 PU2 PU3 N×N 行列 A，B に対し，P 台のプロセッサを使って C = AB を計算 A はブロック行分割，B はブロック列分割されているとする。 結果の C はブロック列分割の形で求めるとする。 各プロセッサは，行列の一部のみを持つ （分散メモリ）。 また，自分の持つ部分のみにアクセスできる。 PU0 PU3 PU0 A0 PU1 A1 × B0 B1 B2 B3 = C0 C1 C2 C3 PU2 A2 PU3 A3

行列乗算の例 (2) 計算の方針 各 CJ を縦方向に P 個のブロック C0J, C1J, …, CP-1,J に分ける。 まず，自分の持つ部分行列のみで計算できる CJJ を計算 その後，他のプロセッサからデータをもらいながら， 他の CIJ を計算 program matmult （配列の確保： AJ,BJ,CJ ,ATMP） （MPIの初期化と環境情報取得。自分のプロセッサ番号をJとする。） （配列AJ,BJの初期化） CJJ = AJ×BJ do I = 1, P-1 （プロセッサ MOD(J-I+P,P)にAJを送る。） （プロセッサ MOD(J+I,P)からAMOD(J+I,P)を受け取り，ATMPに格納。） 　CMOD(J+I,P),J = ATMP×BJ end do （配列CJの出力） （MPIの終了） stop end

3.2　高性能化の技法 基本的な考え方 実行時間の予測 行列乗算の例

基本的な考え方 PU間の負荷分散均等化 PU間通信量の削減 PU間通信回数の削減 キャッシュメモリの有効利用 各PUの処理量が均等になるよう 　　処理を分割 PU間通信量の削減 通信には，データ1個あたり数十サイクルの時間が必要 データ分割の最適化による通信量の削減と通信の隠蔽が重要 PU間通信回数の削減 1回の通信には通常，数百～数千サイクルの起動時間が必要 並列粒度の増大による通信回数の削減が重要 キャッシュメモリの有効利用 キャッシュ PU0 PU1 PU2 PU3 メモリ ネットワーク 並列粒度とは，PU間での同期なしに行える処理の大きさ。

実行時間の予測 演算時間 通信時間 待ち時間 全実行時間 Tcomp = （演算量）／（演算速度） もっとも単純なモデル化 より精密には，キャッシュの影響などを考慮する必要あり 通信時間 Tcomm = （転送回数）×（起動時間） + （転送量）／（転送速度） 待ち時間 Tidle 全実行時間 Texec = Tcomp + Tcomm + Tidle

行列乗算の例 × = アルゴリズム 実行時間の予測 PU0 PU3 PU0 PU1 PU2 PU3 P を増やすと，並列化効率が大きく低下 A，C をブロック行分割，B をブロック列分割とした行列乗算 実行時間の予測 Tcomp = 2N3 / P / s ：P に反比例して減少 Tcomm = (P – 1 ) * (Tsetup + (N / P) * N * 8 / b) = (P – 1 ) * Tsetup + 8 (1 – 1 / P) N2 / b) ：P とともに減少しない Tidle = 0 ：負荷は完全に均等 PU0 PU3 PU0 A0 PU1 A1 × B0 B1 B2 B3 = C0 C1 C2 C3 PU2 A2 PU3 A3 P を増やすと，並列化効率が大きく低下

通信の隠蔽 アイディア 行列 A のデータを，計算に必要な1ステップ前に送ってもらう。 演算とのオーバーラップにより，通信時間を隠蔽できる。 ただし，通信時間 ≦ 演算時間 でないと，隠蔽は不可能 program matmult （配列の確保： AJ,BJ,CJ ,ATMP） （MPIの初期化と環境情報取得。自分のプロセッサ番号をJとする。） （配列AJ,BJの初期化） （プロセッサ MOD(J-1+P,P)にAJを送る。） CJJ = AJ×BJ do I = 1, P-1 （プロセッサ MOD(J+I,P)からAMOD(J+I,P)を受け取り，ATMPに格納。） IF (I < P-1) （プロセッサ MOD(J-I+1+P,P)にAJを送る。） 　CMOD(J+I,P),J = ATMP×BJ end do （配列CJの出力） （MPIの終了） stop end

Fox のアルゴリズム (1) × = × = アルゴリズム A B C プロセッサを2次元に並べ，A，B，C を縦横両方向にブロック分割 プロセッサ (I,J) は，第 K ステップで 　 　　 CIJ += A I, MOD(I+K-1,√P) B MOD(I+K-1,√P), J 　を計算 AII の横方向 ブロードキャスト B0 B1 B2 B3 B0 B1 B2 B3 × = A B C CIJ += A II,B IJ の計算 B0 B1 B2 B3 B0 B1 B2 B3 × =

Fox のアルゴリズム (2) × = × = × = 第2ステップ A B C AI,I+1 の横方向 ブロードキャスト BIJ の縦方向 サイクリック転送 B0 B1 B2 B3 B0 B1 B2 B3 × = A B C CIJ += A I,I+1,B I+1,J の計算 B0 B1 B2 B3 B0 B1 B2 B3 × =

Fox のアルゴリズム (3) アルゴリズム program matmult （配列の確保： AIJ,BIJ,CIJ ,ATMP, BTMP） （MPIの初期化と環境情報取得。自分のプロセッサ番号を(I,J)とする。） （配列AIJ,BIJの初期化） do K = 1, √P 　　IF (J = MOD(I+K,√P)) （配列AIJをATMPにコピー） 　　配列ATMPを行方向のプロセッサ間でブロードキャスト 　　IF (K > 1) THEN 　　　　（プロセッサ (MOD(I-1+P,P),J)にBTMPを送る。） 　　　　（プロセッサ (MOD(I+1,P)からBTMPを受け取る。） 　　ELSE 　　　　（BIJを配列BTMPにコピー） END IF 　CIJ += ATMP×BTMP end do （配列CIJの出力） （MPIの終了） stop end

Fox のアルゴリズムの MPI による実装 × = × = ATMP の横方向ブロードキャスト BTMP の縦方向サイクリック転送 同じ I の値を持つ PU 群に対してコミュニケータ（PUのグループ）を定義 コミュニケータを用い，MPI_BCAST によりブロードキャスト BTMP の縦方向サイクリック転送 MPI_SEND を用いて１つ上の PU にデータを送信 MPI_RECV を用いて1つ下の PU からデータを受信 × = B0 B1 B2 B3 B0 B1 B2 B3 × =

Fox のアルゴリズムの性能モデル 実行時間の予測 P を増やすと，通信時間も 1/√P で減少 Tcomp = 2N3 / P ：P に反比例して減少 Tcomm = √P * log2(√P) * (Tsetup + N2 / P * 8 / b) ：ATMP の転送 + (√P – 1 ) * (Tsetup + N2 / P * 8 / b) ：BTMP の転送 Tidle = 0 ：負荷は完全に均等 P を増やすと，通信時間も 1/√P で減少 P を増やしたときの並列化効率の低下が小さい。

より効率的な行列乗算 通信の隠蔽 より効率的なアルゴリズムの利用 ATMP の転送を，計算の1ステップ前に行う。 通信用バッファがもう1個必要 ノンブロッキング通信（MPI_ISEND，MPI_IRECV）を使うと，より効果的 より効率的なアルゴリズムの利用 SUMMA (Scalable Universal Matrix Multiplication) LAPACK Working Notes 96 http://www.netlib.org/lapack/lawns/downloads/

連立一次方程式の高性能解法 （密行列の場合） 連立一次方程式の高性能解法　（密行列の場合） １.　LU分解 ２.　LU分解のブロック化 ３.　LU分解の並列化 本発表では，はじめに，研究の背景を述べてから，スパースソルバの概要，並列化手法，そして本研究で工夫した点の一つであるRISCプロセッサ向けの最適化についてご説明します。最後に，並列計算機SR2201上での性能評価とまとめを述べます。

ガウスの消去法の性能 n=1000のときの性能は250MFLOPS程度