線形方程式求解アルゴリズムに対する求解性能データに基づいた体系的な性能比較– 行列・固有値問題の解法とその応用 –-第２回研究会(2006.11.16)

伊藤祥司(筑波大学 大学院システム情報工学研究科)

発表内容
1. イントロダクション- モチベーション- 本研究の目的

2. 本研究のアプローチ- 性能評価の情報システム- システム環境，今回の評価方法- Survey chart（Visual Evaluation）

3. 実験結果，考察- Solver ? orPreconditioning ?
- I+S 前処理の効果について

4. まとめはじめに
自然現象や工学現象の解析では，数値シミュレーションが行われることもしばしばであり，それらを記述する方程式は，多くの場合，大規模な線形（連立１次）方程式

Ax = b

A : n×n 係数行列，x, b : n次ベクトルを解くことに帰着される．
シミュレーションに要する計算時間の大半が，この計算に費やされるため，速く，正確に解くことが重要である．

〔モチベーション〕どのアルゴリズムを選択すれば良いのだろうか？
直接解法
ガウス消去法に基づく解法,LU , Cholesky
others

アルゴリズムを体系的に性能比較した情報が，十分に提供されていない！
反復解法
非定常反復法(Krylov部分空間法)

更には :解法と併用する前処理も !
ILU, Point Jacobi,Hybrid, I+S, SA-AMG
othersConjugate Gradient
BiCG , CGS
Bi-CGSTAB
BiCGstab(l )GMRES , GPBiCG
others

定常反復法
Jacobi
Gauss-SeidelSOR
othersa Survey and Evaluation System
on Numerical Algorithms

本研究の目的
数値計算アルゴリズムの体系的な性能評価と特性分析
求解アルゴリズムが生成する様々なデータや解くべき問題の性質を示すデータその他の注目すべきファクター等に対して，品質管理，統計分析，データマイニングや可視化評価などを用いて，実際の特性や傾向などを分析する．

評価方法の一例（イメージ）
前処理 A

行列の何らかの特性

残差ノルム
反復回数+++++

解法 A ++
前処理 B +++

アルゴリズムからの情報

アルゴリズムの個々の性能評価:

アルゴリズムの体系的な性能評価:

１つの問題（係数行列） に対して，複数のアルゴリズムの収束の様子を評価してきた．

様々な問題の，ある特性と，アルゴリズムとの間にはどのような関係が存在するだろうか？

発表内容
1. イントロダクション- モチベーション- 本研究の目的

2. 本研究のアプローチ- 性能評価の情報システム- システム環境，今回の評価方法- Survey chart（Visual Evaluation）

3. 実験結果，考察- Solver ? orPreconditioning ?
- I+S 前処理の効果について

4. まとめ情報システム（概念図）
Matrices
行列 A

Algorithms

行列 B

Bi-CGSTAB

行列 C
その他
利用者はWWW経由で閲覧可能

計算サーバ

BiCG

計算により得られたデータ

BiCGstab(l )GMRESothersDB

ﾃﾞｰﾀを解析

データ解析
サーバ

Webサーバ
(DBへのアクセス /ﾃﾞｰﾀの可視化)

システム環境(計算サーバ)- Matrices (Test problems):
線形方程式向きの52 種類の行列 （ Matrix Marketの行列 ）各々のRHSは b=Ax, x=(1.0 , … , 1.0)を用いて生成した- Algorithms (Numerical Calculation Library):Lis [Lis-AMG-1.0.1 seq. ver.] (小武守，西田，長谷川：CREST SSI prj.)
12 解法 × 8 前処理(“前処理無し”も含める )MachineSun Fire V880

CPUUltraSPARC III (900 MHz)

Memory size

16 GB

Operating System

Solaris 9Compiler

Sun Workshop 6 (cc, f90)

Outline of Lis libraryLis (a Library of Iterative Solvers for linear systems)- 線形方程式求解アルゴリズムの逐次版，並列版（MPI, OpenMP）.
- １１種類のデータ形式をサポート. (CRS, CCS, MSR, DIA, …)Component of Lis: [Lis-AMG-1.0.1 modified ver. (ver.1.0.2とほぼ同じ) ](† Lis の特徴：国内の研究成果)
反 復 解 法
非定常反復解法
CG, BiCG, CGS, BiCGStab, BiCGStab(l ),GPBiCG†, TFQMR, Orthomin(m), GMRES(m)定常反復解法Jacobi, Gauss-Seidel, SOR

前 処 理
非定常反復解法向きNone, (Point) Jacobi, ILU(k), SSOR, Hybrid†,
I+S†, SAINV, SA-AMG†
定常反復解法向き
None, I+S†

評価条件
- 最大反復回数:
行列のサイズを設定- 反復法の初期ベクトル:x0 = 0- 収束判定条件:

rk 2 b 2 ≤ 1.0 ×10−12- ベクトル r ： アルゴリズム中の残差ベクトル ,rk +1 = rk − α k Apk（CG法の例）- 数値解が有効かどうかの判定 ： 真の残差ベクトルrˆ = b − Axˆ ,x̂ : 数値解rˆ 2 b 2 ≤ 1.0 ×10−8

Survey chart

Matrices(Grouping by solver)

Algorithms (Combination of Solver and Preconditioning)チャート内の配色各行列から生成される各々の線形方程式に対して：最速のアルゴリズム
Score =他のアルゴリズム* 10(収束までの所要反復回数やCPU時間について)“*” “真の残差の判定では収束していない”(ただし，アルゴリズム中の残差では収束している)“.” “最大反復回数に到達しても収束しない，または，ブレークダウンで反復が停止した”(要するに，アルゴリズム中の残差でも収束していない)発表内容
1. イントロダクション- モチベーション- 本研究の目的

2. 本研究のアプローチ- 性能評価の情報システム- システム環境，今回の評価方法- Survey chart（Visual Evaluation）

3. 実験結果，考察- Solver ? orPreconditioning ?
- I+S 前処理の効果について

4. まとめSurvey chart on Iter. Num.
(Grouping by Solver)

CG

BiCG

CGS

BiCGStab

BiCGStab(l=2)

GPBiCGTF
QMR

Orthomin

GMRES

Jacobi

GaussSeidel

SOR

Survey chart on Iter. Num.
(Grouping by Preconditioning)

None

Point
JacobiHybridILU(k)(IC(k))

SSOR

SAINV

I+SSAAMG

考察（その１）収束性に対する効果は，解法の違いよりも前処理の違いによる影響の方が大きいことが確認できる（本実験の条件下で）．

For NonSymmetric systems
CPUtime
None

Point
JacobiILU(k)

Iter.
Num.

SSOR

HybridI+SSAINVSAAMG

For Symmetric systems
CPUtime
None

Point
JacobiILU(k)

SSOR

Iter.
Num.

SAINV

Hybrid(IC(k))

I+SSAAMG

I+S 前処理
I+S 前処理 ： Jacobi法 Gauss-Seidel法向けに開発された
A.D.Gunawardena, S.K.Jain and L.Snyder (1991)

行列 S :

係数行列 :⎡ a11⎢a⎢ 21
A=⎢ M
⎢⎢ M
⎢⎣an1

a12

a13 L

a22
M
M
an 2a23
O
M
an 3

L
LO
L

a1n ⎤a2 n ⎥⎥M ⎥⎥an −1n ⎥ann ⎥⎦

Ax = b(m)(m)( I + S ) Ax = ( I + S )b

⎡0 − a12⎢00⎢
S = ⎢00⎢M⎢M⎢⎣000

L− a23 L0OM0OL⎤
M ⎥⎥M ⎥⎥− an −1n ⎥0 ⎥⎦0m 本の副対角要素

岡山理科大： 仁木，河野，小武守，他

I : 単位行列・I+S 前処理に対する数々の改良

m=1: original, m=3: Lis default

・Krylov部分空間法に対する適用考察（その２）対称行列に対する I+S 前処理の効果：
一般に，行列 A が対称であっても，前処理 I+S を施した系は非対称である.
通常，係数行列 A が対称であるとき，Krylov部分空間法の中で共役勾配(CG)法がより良い（対称性を保持）と考えられている．

しかし，本研究の結果からは，対称系に対しても，CG法に基づくアルゴリズムより，I+S前処理とCG法以外の解法の組合せの方が速く収束する事例も確認できる．

まとめ1) 線形方程式の求解アルゴリズムに対する体系的な性能評価を行うためのシステムを提案した．
性能情報検索ｼｽﾃﾑ (暫定URL: ﾃﾞｰﾀ，機能の一部を公開)
http://mma.cs.tsukuba.ac.jp/~itosho/sesna/2) Lisを用いた体系的な性能比較の結果から，反復解法の収束に対する影響は，解法よりも前処理の方が大きいことが確認された．
3) 特に，I+S 前処理は，非対称系だけでなく対称系に対しても効果があることが確認された．