Transcript IDR(S)?? - 筑波大学大学院 システム情報工学研究科 コンピュータ

一般化Bi-CGSTAB(s, L)
(=一般化IDR(s, L))
東京大学大学院 情報理工学系研究科
数理情報学専攻
谷尾 真明 杉原 正顯
1
発表の流れ


概要
IDR(s)法, GIDR(s,L)法の導出





IDR(s)法
GIDR(s,L)法
GBi-CGSTAB(s,L)法の導出

(応用数理年会)
(RIMS研究集会)
GBi-CG(s)法
GBi-CGSTAB(s,L)法(=GIDR(s,L)法)
まとめ
2
発表の流れ


概要
IDR(s)法, GIDR(s,L)法の導出





IDR(s)法
GIDR(s,L)法
GBi-CGSTAB(s,L)法の導出

(応用数理年会)
(RIMS研究集会)
GBi-CG(s)法
GBi-CGSTAB(s,L)法(=GIDR(s,L)法)
まとめ
3
Krylov部分空間法
概要
Krylov部分空間法
Bi-CG
CGS
Bi-CGSTAB
Bi-CGSTAB(L)
GP-BiCG
Lanczos原理
(短い漸化式)
GMRES
FOM
GCR
Orthomin
Arnoldi原理
(長い漸化式)
IDR(induced dimension reduction)(s)[1](2007)
(短い漸化式)
[1]P. Sonneveld and M. B. van Gijzen, IDR(s): A Family of Simple and Fast Algorithms
for Solving Large Nonsymmetric Systems of Linear Equations, Reports Depart. Appl.
4
Math. Anal., REPORT 07-07, Delft Univ., Tech., 2007.
Lanczos原理に基づくKrylov部分空間法
vs IDR(s)法
概要
一番重たい計算＝行列とベクトルの積
•
Krylov部分空間法
高々2N回行えば収束(理論)
•
IDR(s)法
A
u
v
高々N+N/s回行えば収束(理論)
行列ベクトル積
s>1の時、既存の方法より少ない
実際の数値計算においても
早く収束する！
１回
5
IDR(s)法の一般化
IDR(s)法
計算量
加速多項式
の次数
概要
GIDR(s,L)法[2]
N+N/s
1
N+N/s
L
歪対称に近い係数行列と
相性が良くない
[2]谷尾真明 杉原正顯, GIDR(s,L): 一般化IDR(s),
応用数理学会2008年度年会, pp. 411—412, 千葉.
6
既存のKrylov部分空間法と
IDR(s), GIDR(s,L)の関係図


概要
IDR(1) = Bi-CGSTAB
GIDR(1,L) = Bi-CGSTAB(L)
加速多項式の次数
0
1
s
1
L
Bi-CG
Bi-CGSTAB
=IDR(1)
Bi-CGSTAB(L)
=GIDR(1,L)
？
IDR(s)
GIDR(s,L)
拡張
拡張
7
既存研究
概要
Bi-CG(s)⇒Bi-CGSTAB(s)(=IDR(s))
(Sleijpen, 2008)

加速多項式の次数
0
1
s
Bi-CG
Bi-CG(s)
1
L
Bi-CGSTAB
=IDR(1)
Bi-CGSTAB(L)
=GIDR(1,L)
IDR(s)
=Bi-CGSTAB(s)
GIDR(s,L)
拡張
8
本研究

概要
GBi-CG(s)⇒GBi-CGSTAB(s,1)(=IDR(s))
⇒GBi-CGSTAB(s,L)(=GIDR(s,L))
加速多項式の次数
0 拡張
1
s
Bi-CG
Bi-CG(s)
GBi-CG(s)
1
L
Bi-CGSTAB
=IDR(1)
Bi-CGSTAB(L)
=GIDR(1,L)
IDR(s)
=Bi-CGSTAB(s)
=GBi-CGSTAB(s,1)
GIDR(s,L)
GBi-CGSTAB(s,L)
9
発表の流れ


概要
IDR(s)法, GIDR(s,L)法の導出





IDR(s)法
GIDR(s,L)法
GBi-CGSTAB(s,L)法の導出

(応用数理年会)
(RIMS研究集会)
GBi-CG(s)法
GBi-CGSTAB(s,L)法(=GIDR(s,L)法)
まとめ
10
IDR(s)法
IDR定理(1/2)

定義:
: span{
}
と直交する空間
: N次元ベクトル空間
(注意)
は非ゼロのスカラー量
11
IDR定理(2/2)

IDR(s)法
定理
(genericな条件の下で)
-s次元
-s次元
・
12
IDR(s)法
IDR(s)法の概要

残差ベクトルが入れ子になっている空間列
の中を
と潜っていくように更新.
ri  b  Axi
r1
残
差
ベ
ク
ト
ル
r2
r0 r3
N+N/sMATVECsで0に到達
13
発表の流れ


概要
IDR(s)法, GIDR(s,L)法の導出





IDR(s)法
GIDR(s,L)法
GBi-CGSTAB(s,L)法の導出

(応用数理年会)
(RIMS)
GBi-CG(s)法
GBi-CGSTAB(s,L)法(=GIDR(s,L)法)
まとめ
14
IDR(s)法の不満点


GIDR(s,L)法
係数行列が歪対称行列に近いと収束性
悪化.
IDR定理自身が持つ問題点.
加速多項式に相当
→高次にしたい！
15
GIDR(s,L)法
数値実験による確認

IDR(3) (1次の加速多項式)
vs Bi-CGSTAB(3) (3次の加速多項式)
3次元対流問題
uxx  u yy  uzz 1000ux  F
離散化
Ax  b
・125000×125000
・歪対称行列に近い
Bi-CGSTAB(3)
IDR(3)
16
GIDR(s,L)法
一般化IDR定理(1/2)
一般化

定義:
: span{
}
: N次元ベクトル空間
高次
:
17
GIDR(s,L)法
一般化IDR定理(2/2)

定理
(genericな条件の下で)
-sL次元 -sL次元
・
18
GIDR(s,L)法
GIDR(s,L)のアルゴリズム
1. repeat until
2. For i = 0,1,…,L-1
3.
For j = 1,…,s
4.
Solve
5.
for
(k=0,…,i)
6.
7.
8.
9.
10.
11.
12.
13.
14.
end
15. Select
16. For i=1,…,L
17.
18.
19.
Compute
end
Solve
for
20. end
21.
22.
23. end repeat
残差の更新の際に
補助ベクトルを導入している
19
GIDR(s,L)法
数値実験による確認

IDR(3) (1次の加速多項式)
vs GIDR(3,3) (3次の加速多項式)
3次元対流問題
uxx  u yy  uzz 1000ux  F
離散化
Ax  b
・125000×125000
・歪対称行列に近い
GIDR(3,3)
IDR(3)
20
発表の流れ


概要
IDR(s)法, GIDR(s,L)法の導出





IDR(s)法
GIDR(s,L)法
GBi-CGSTAB(s,L)法の導出

(応用数理年会)
(RIMS研究集会)
GBi-CG(s)法
GBi-CGSTAB(s,L)法(=GIDR(s,L)法)
まとめ
21
関係図
M
= : Mathematically equivalent
加速多項式の次数
1
s
0
1
L
Bi-CG
Bi-CGSTAB
Bi-CGSTAB(L)
IDR(s)
GIDR(s,L)
＝Bi-CG(s)
M
=GBi-CG(s)
(ML(s)BiCG)
M
=IDR(1)
M
=Bi-CGSTAB(s)
M
=GBi-CGSTAB(s,1)
(ML(s)BiCGSTAB)
M
=GIDR(1,L)
GBi-CGSTAB(s,L)
22
GBi-CG(s)
Bi-CG法

Krylov 部分空間

反復
条件
shadow residual
23
GBi-CG(s)
Bi-CG法の1反復
s. t.
s. t.
残差が満たす性質
24
Bi-CG法の一般化
shadow residualの高次化
GBi-CG(s)
Bi-CG
shadow
residual
GBi-CG(s)
N
・・・
1
s
条件
span
N
Bi-CG法の一般化
shadow residualの高次化



GBi-CG(s)
ML(s)BiCG(Yeung and Chan, 1999)
Bi-CG(s) (Sleijpen, 2008)
GBi-CG(s)
加速多項式の次数
shadow residual
0
1
L
1
Bi-CG
Bi-CGSTAB
Bi-CGSTAB(L)
s
Bi-CG(s)
GBi-CG(s)
ML(s)BiCG
26
GBi-CG(s)
GBi-CG(s)の1反復
s. t.
For j = 1,…,s
set
s. t.
end
残差が満たす性質
27
GBi-CG(s)法の計算量

1反復あたりの行列ベクトル積の演算



GBi-CG(s)
補助ベクトル ・・・・・・・・・s
shadow residual ・・・・・・s
×s
収束するのに必要な反復回数
(genericな条件下で)
total 2s × N/s = 2N
MATVECs
(Bi-CG法, Bi-CGSTAB法などと同じ)
28
発表の流れ


概要
IDR(s)法, GIDR(s,L)法の導出





IDR(s)法
GIDR(s,L)法
GBi-CGSTAB(s,L)法の導出

(応用数理年会)
(RIMS研究集会)
GBi-CG(s)法
GBi-CGSTAB(s,L)法(=GIDR(s,L)法)
まとめ
29
GBi-CGSTAB(s,L)
拡張のアイデア

Bi-CG⇒Bi-CGSTAB(L)と同じ構造を利用すること
によりGBi-CG(s)⇒GBi-CGSTAB(s,L)を達成！
加速多項式の次数
1
s
0
1
L
Bi-CG
Bi-CGSTAB
Bi-CGSTAB(L)
IDR(s)
GIDR(s,L)
＝Bi-CG(s)
M
=GBi-CG(s)
(ML(s)BiCG)
M
=IDR(1)
M
=Bi-CGSTAB(s)
M
=GBi-CGSTAB(s,1)
(ML(s)BiCGSTAB)
M
=GIDR(1,L)
GBi-CGSTAB(s,L)
30
加速多項式の仕組み
(Bi-CG) (1/3)

GBi-CGSTAB(s,L)
Bi-CG⇒Bi-CGSTAB
⇒Bi-CGSTAB(L)
加速多項式の次数
1
s
0
1
L
Bi-CG
Bi-CGSTAB
Bi-CGSTAB(L)
IDR(s)
GIDR(s,L)
＝Bi-CG(s)
M
=GBi-CG(s)
(ML(s)BiCG)
M
=IDR(1)
M
=Bi-CGSTAB(s)
M
=GBi-CGSTAB(s,1)
(ML(s)BiCGSTAB)
M
=GIDR(1,L)
GBi-CGSTAB(s,L)
31
加速多項式の仕組み
(Bi-CG) (2/3)
GBi-CGSTAB(s,L)
・Bi-CG法における任意性
s. t.
s. t.
32
加速多項式の仕組み
(Bi-CG) (3/3)
Bi-CG
Bi-CG
+
加速多項式
残差
加速多項式
GBi-CGSTAB(s,L)
shadow residual
(fixed)
33
加速多項式の仕組み
(GBi-CG(s)) (1/3)

GBi-CGSTAB(s,L)
M
GBi-CG(s)⇒GBi-CGSTAB(s,1)(=IDR(s))
⇒GBi-CGSTAB(s,L)(=GIDR(s,L))
加速多項式の次数
1
s
0
1
L
Bi-CG
Bi-CGSTAB
Bi-CGSTAB(L)
IDR(s)
GIDR(s,L)
＝Bi-CG(s)
M
=GBi-CG(s)
(ML(s)BiCG)
M
=IDR(1)
M
=Bi-CGSTAB(s)
M
=GBi-CGSTAB(s,1)
(ML(s)BiCGSTAB)
M
=GIDR(1,L)
GBi-CGSTAB(s,L)
34
加速多項式の仕組み
(GBi-CG(s)) (2/3)
GBi-CGSTAB(s,L)
・GBi-CG(s)における任意性
For i = 1,…,s
set
s. t.
end
残差が満たす性質
35
加速多項式の仕組み
(GBi-CG(s)) (3/3)
残差
GBi-CGSTAB(s,L)
shadow residual
GBi-CG(s)
GBi-CG(s)
+
加速多項式
加速多項式
(fixed)
36
GBi-CGSTAB(s,L)
GBi-CGSTAB(s,L)=GIDR(s,L)
加速多項式の次数
1
s
0
1
L
Bi-CG
Bi-CGSTAB
Bi-CGSTAB(L)
IDR(s)
GIDR(s,L)
＝Bi-CG(s)
M
=GBi-CG(s)
(ML(s)BiCG)
M
=IDR(1)
M
=Bi-CGSTAB(s)
M
=GBi-CGSTAB(s,1)
(ML(s)BiCGSTAB)
M
=GIDR(1,L)
GBi-CGSTAB(s,L)
37
GBi-CGSTAB(s,L)
GBi-CGSTAB(s,L)
GIDR(s,L)のアルゴリズム
1. repeat until
2. For i = 0,1,…,L-1
3.
For j = 1,…,s
4.
Solve
5.
for
(k=0,…,i)
6.
7.
8.
9.
10.
11.
12.
13.
再掲
14.
end
15. Select
16. For i=1,…,L
17.
18.
19.
Compute
end
Solve
for
20. end
21.
22.
23. end repeat
補助ベクトル
⇔
GBi-CG(s)における補助ベクトル
38
GBi-CGSTAB(s,L)
GBi-CGSTAB(s,L)の計算量
k反復後の残差
k反復後の
shadow residual
GBi-CG(s)
計算量
2N
Matvec
ks
Matvec
k = [N/s]で収束
GBi-CG(s)
+
加速多項式
計算量
N+N/s
ks
(fixed)
Matvec
ks + k
Matvec
0
39
発表の流れ


概要
IDR(s)法, GIDR(s,L)法の導出





IDR(s)法
GIDR(s,L)法
GBi-CGSTAB(s,L)法の導出

(応用数理年会)
(RIMS研究集会)
GBi-CG(s)法
GBi-CGSTAB(s,L)法(=GIDR(s,L)法)
まとめ
40
結論

GIDR(s,L)法を, Bi-CG法の拡張として解釈
出来た(GBi-CGSTAB(s,L)法)
1
s
0
1
L
Bi-CG
Bi-CGSTAB
Bi-CGSTAB(L)
IDR(s)
GIDR(s,L)
＝Bi-CG(s)
M
=GBi-CG(s)
(ML(s)BiCG)
M
=IDR(1)
M
=Bi-CGSTAB(s)
M
=GBi-CGSTAB(s,1)
(ML(s)BiCGSTAB)
M
=GIDR(1,L)
GBi-CGSTAB(s,L)
41
今後の課題

前処理を含めたGBi-CGSTAB(s, L)法の提案.

パラメータsとLの決定法.
0
shadow residual
1
s
加速多項式の次数
1
L
Bi-CG
Bi-CGSTAB
Bi-CGSTAB(L)
Bi-CG(s)
GBi-CG(s)
ML(s)BiCG
IDR(s)
Bi-CGSTAB(s)
GBi-CGSTAB(s,1)
ML(s)BiCGSTAB
GIDR(s,L)
GBi-CGSTAB(s,L)
42