応用数理工学特論　第5回 計算理工学専攻　張研究室 山本有作

前回の概要 連立一次方程式の高性能解法 （密行列の場合） LU分解 LU分解のブロック化 LU分解の並列化（共有メモリ型並列計算機） 連立一次方程式の高性能解法　（密行列の場合） LU分解 LU分解のブロック化 LU分解の並列化（共有メモリ型並列計算機） LU分解の並列化（分散メモリ型並列計算機）

今回の概要 連立一次方程式の高性能解法 （疎行列の場合） 差分法により生じる疎行列 疎行列向け連立一次方程式解法の分類 帯ガウス法 連立一次方程式の高性能解法　（疎行列の場合） 差分法により生じる疎行列 疎行列向け連立一次方程式解法の分類 帯ガウス法 スパースソルバ

差分法により生じる疎行列 ポアソン方程式 応用 物理的意味 熱伝導 静電場 流体力学 応力とひずみの解析 ある点でのu の値 = その周囲でのu の値の平均値 + 外力の項

大規模計算の例 ３次元の応力・ひずみ解析

簡単な例 ２次元正方領域での熱伝導方程式 y 1 u(0, y) = 0 u(1, y) = 0 u(x, 0) = 0 1

２次元ポアソン方程式の離散化 ５×５の２次元格子 差分法で離散化 → 疎行列 格子 行列 1 2 3 4 5 6 7 8 9 10 11 差分法で離散化　→　疎行列 格子 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 行列

疎行列向け連立一次方程式解法の分類 直接法 --- ガウスの消去法に基づく方法 反復法 --- 解を逐次的に改良する方法 帯ガウス法 スパースソルバ 反復法 --- 解を逐次的に改良する方法 定常反復法 ヤコビ法，ガウス=ザイデル法，SOR法 クリロフ部分空間法 CG法（共役勾配法） CGS法，Bi-CGSTAB法，GPBi-CG法 GMRES法

直接法の特徴 長所 短所 安定性 必要なメモリ量・計算時間が事前に予測可能 行列が対称正定値であれば，必ず解が求まる。 必要なメモリ量・計算時間が事前に予測可能 係数行列が同じで右辺ベクトルのみ異なる複数の連立一次方程式を解くのに有利 短所 必要なメモリ量が大きい フィルインのため 反復法がうまくいく場合は，計算時間の面で不利

反復法の特徴 長所 短所 必要なメモリ量が小さい。 うまくいく場合は非常に高速 収束しない場合や，収束に非常に時間のかかる場合がある。 元の行列（+ 数本の作業用ベクトル）のみを使って計算可能 うまくいく場合は非常に高速 短所 収束しない場合や，収束に非常に時間のかかる場合がある。

反復法の収束特性の例 ヘルムホルツ方程式 張研で開発 されたGPBi-CG法→

帯ガウス法 帯行列に対してガウス消去法を適用すると，非ゼロ要素は帯の内部にしか現れないことに着目 帯の内部のみを記憶・演算して演算量とメモリ量を節約

２次元ポアソン方程式のオーダリング ５×５の２次元格子 分割を２段階行い，領域を４個の部分領域に分割 格子 行列 ND法による オーダリング 1 2 3 4 5 ND法による オーダリング 1 3 9 5 7 第１段階での セパレータ 6 7 8 9 10 2 4 10 6 8 21 22 23 24 25 11 12 13 14 15 第２段階での セパレータ 16 17 18 19 20 11 13 19 15 17 21 22 23 24 25 12 14 20 16 18 行列 ND法による オーダリング 各対角ブロック の分解は独立に 実行可能

一般の場合のオーダリング Nested Dissection 法によるオーダリング （有限要素法から生じる疎行列の場合） セパレータ Nested Dissection 法によるオーダリング 　 （有限要素法から生じる疎行列の場合） セパレータと呼ばれる節点集合を取り除くことにより，メッシュを２つの部分に分割 この分割に従って行列の行・列を並び替え，行列を縁付きブロック対角形に変形 この処理を各対角ブロックに対して再帰的に繰り返す オーダリングによる並列性抽出 コレスキー分解のステップでは，各対角ブロックを並列に分解可能 各対角ブロックの演算負荷が均等で，かつ縁の部分が狭い（セパレータに属する節点の個数が小さい）分割が望ましい 有限要素法メッシュ 置換後の行列 1 2 3 4 5 6 7 7 ３段階の分割後の行列

消去木と演算の並列性 消去木とは ND法でオーダリングした行列に対する消去木 疎行列のコレスキー分解における消去演算の依存関係を表現するグラフ 消去木の各節点は行列 L の各列に対応 ND法でオーダリングした行列に対する消去木 消去木と行列の対応関係 鎖の部分 ⇔ 縁の部分 部分木 ⇔ 対角ブロック 共通部分を持たない２つの部分木に対する分解演算は並列に実行可能 再帰的縁付きブロック 対角行列 対応する消去木

消去木を用いたデータ割り当て subtree-to-subcube 割り当て 消去木が不均等な場合 鎖の部分では，周期的に節点をプロセッサに割り当てる 分岐点でプロセッサ群を２分割し，各部分木をそれぞれのプロセッサ群に割り当てる。この操作を再帰的に行う 消去木が不均等な場合 単純な subtree-to-subcube 割り当てでは大きな負荷不均衡が生じる場合あり subtree-to-subcube 割り当ての改良が必要

コレスキー分解 処理 アルゴリズム 各プロセッサに分散された行列 A を入力として，並列にコレスキー分解を行う 対角ブロックの分解と縁の部分の計算を，小さいレベルのブロックから順に再帰的に繰り返す 入力行列 消去木 ４プロセッサで協調して消去演算 ２プロセッサで協調して消去演算 各プロセッサで独立に消去演算（通信なし）

スパースソルバの処理フロー START オーダリング シンボリック分解I 各プロセッサへの データ割り当て決定 データ分散 置換行列Pの作成。演算量削減と並列性抽出が目標 分解後の非零要素決定，演算量・所要メモリ量算定， 消去木作成 シンボリック分解I 各プロセッサへの データ割り当て決定 消去木に基づき，負荷均等化・通信最小化を達成する ように行列の列のプロセッサへの割り当てを決定 データ分散 上記の割り当て方法に基づき，行列データを分散 シンボリック分解II 分解後の非零要素を格納するデータ構造を作成 コレスキー分解 前進後退代入 データ収集 END

試作版スパースソルバの性能評価結果 概要 計算機環境 以上のアルゴリズムに基づく並列プログラムを開発し，２種のPCクラスタ上で性能評価を行った 計算機環境 Alphaクラスタ（mellon） Alpha EV5 533Mhz ×8ノード RAM 128MB／ノード 100base-tx 接続 コンパイラは g77 –O3 を使用 Power PC G5クラスタ（muscat） PowerPC G5 2.0Ghz × 8ノード RAM 2GB／ノード Gigabit Ether接続

例題 Harwell Boeing Matrix Collectionの５種の行列と，３次元ポアソン方程式を離散化して得られる行列に対し，コレスキー分解の性能を評価 Title Discipline Size(N) Nonzero 演算量(Gflop) 3dp24 ３次元ポアソン 243 41472 1.15 3dp32 323 98304 6.59 BCSSTK17 Elevated pressure vessel 10974 428650 0.24 BCSSTK18 R.E. Ginna Nuclear Power Station 11948 149090 0.11 BCSSTK29 Buckling model of a Boeing 767 rear pressure bulkhead 13992 619488 0.66 BCSSTK30 Statics module of an off-shore generator platform 28924 1036208 1.37 BCSSTK31 Statics module of an automobile component 35538 608502 1.87

１プロセッサでの実行時間（コレスキー分解部） Alphaプロセッサでの実行時間 G5プロセッサでの実行時間 G5 は Alpha の10倍程度高速

並列実行時の加速率 3dp24，32，bcsstk30，31は比較的よい結果 bcsstk29は加速率が低い

消去木による分析 消去木を用いて，各レベルでの負荷分散と演算・通信の時間比を調べ，性能評価結果を分析 3dp24の場合 負荷分散は比較的良い 　　→　加速率良好 上のレベルでの演算が多く， 　 独立に実行可能な演算が 　 少ない 　　→　これを改善すれば，並列 　　　　 化効率を更に改善可能 レベル0： 全PUで協調 レベル1： 4PUで協調 レベル2： 2PUで協調 レベル3： 各PUが独立に計算 演算時間 通信時間 面積が時間の長さを表す

消去木による分析（続き） bcsstk30 bcsstk29 消去木を用いて並列性能を解析し，性能改善に役立てることが可能 ・ 独立に演算できる部分が多い 　　→　加速率良好 ・ 負荷分散はやや不均衡 　　→　これを改善すれば，並列化 　　　　 効率を更に向上可能 ・ 演算のほとんどは最上位レベルで， 　独立に演算できる部分が少ない 　　→　並列化効率の改善には，より 　　　　 良いオーダリングが必要 消去木を用いて並列性能を解析し，性能改善に役立てることが可能

参考文献 山本有作： “疎行列連立一次方程式の直接解法”，計算工学，Vol. 11, No. 4, pp. 14-18 (2006). 　　http://www.na.cse.nagoya-u.ac.jp/~yamamoto/work.html

グループ発表について 概要 スケジュール 最後の２回の授業を使って行う。 ２～３人のグループで１つの発表を行う。 講義で勉強したことをもとに，自分たちで数値実験を行った結果を発表。あるいは論文を読んでまとめて発表。 スケジュール グループ決定：　6/19 まで 発表テーマ決定：　6/26 まで 発表：　7/10，7/17 （各グループ15～20分程度）

発表テーマの例 自分の取り組んでいるシミュレーションのプログラムを，高性能化あるいは並列化した結果の報告 講義で学んだ行列計算アルゴリズムの性能評価 分散メモリ型計算機上でのLU分解の並列化 スパースソルバの性能評価 固有値計算アルゴリズムの性能評価，など 高性能計算に関する論文を読み，まとめて報告 その他