|TOP Page|


Algebraic Multigrid Methods

AMG法、代数的マルチグリッド法、代数的多重格子法




Last Update:2009年05月15日


目次

概要

AMG法はMultigrid法の一種の方法であり、直交格子だけでなく、非構造格子や一般的な疎行列も解けるように幾何マルチグリッド法を拡張したものである。マルチグリッド法と同じく、多くの問題では反復数がメッシュサイズに依存しないという特徴があり、現在もっとも高速な連立一次方程式の解法の1つといえる。

基本的な手続き

AMGが幾何的マルチグリッド法と違う点は、

である。AMG法では行列を解く準備として、

  1. Coarse Grid/Fine Gridを決める(F/C Splitting)
  2. Prolongation行列Pを作る
  3. Restriction行列R,CoarseGridのオペレータA_cを作る

ProlongagionとRistrictionとCoarseGridオペレータ

Smootherで収束が遅い誤差の性質

Coarse Grid CorrectionはSmootherの収束が遅いような誤差を落とすことが目的である。この誤差を取り除くためにはProlongationオペレータの値域にこの誤差が含まれていなければならない。よってこのSmootherの収束が遅いような誤差がどのような性質を持っているのかを調べてそのような誤差を表現できるようにProlongationのオペレータを作成する。Smootherの収束が遅いような誤差は次のような性質を持っている。

以下この上の性質について詳しく述べる

残差が小さい誤差

古典的反復法における残差は次のように更新される

\b{e}_{k+1}=(\b{I}-\b{Q}^{-1}\b{A})\b{e}_k

古典的反復法で収束が遅い場合は\b{e}_{k+1}\simeq\b{e}_kであるから、

(\b{I}-\b{Q}^{-1}\b{A})\b{e}_k\simeq\b{e}_k

\b{e}_k-\b{Q}^{-1}\b{r}_k\simeq\b{e}_k

\b{Q}^{-1}\b{r}_k\simeq0

つまり、\b{r}_k\simeq0

Smootherで収束が遅い誤差は残差が小さいということがいえる。

CoarseGridCorrectionではこのようなSmootherの効率がわるい、残差が小さいという性質をもった誤差を効率よく小さくすることを目標としている。そのため、残差が小さいという性質は非常に重要で、この単純な性質をもった誤差が値域にはいるようにProlongationのオペレータを決め、CoarseGridCorrectionを行う。

このような小さい残差を生むような誤差は代数的に滑らかな誤差(Algebraically Smooth error)という。問題によって代数的に滑らかな誤差でも幾何的には振動している可能性がある。またこのように小さな残差を生むような誤差の別の呼び方として値が0に近いことからNear Kernel errorともいう。

残差が小さい誤差(もう少し正確に)

上の議論は大雑把で、小さい残差とはいっても何と比較して小さいのかわからない。もっと正確には、2つの誤差ノルムに次の関係が成り立つとき残差が小さいという。

(Ae,D^{-1}Ae)<<(e,Ae)つまり、(r,D^{-1}r)<<(e,r)

成分で書くと

\sum_{i=1}^n\frac{r_i^2}{a_{ii}}<<\sum_{j=1}^ne_jr_j

上の式は和で表されているが、平均的な性質として

|r_i|<<a_{ii}|e_j|

となることが期待される。

強い接続に沿って誤差の変化が小さいという性質

さて、残差が小さいつまり、(Ae,D^{-1}Ae)<<(e,Ae)であるような誤差はどのような性質を持つのかを調べてみよう。

次の関係が成り立つ。

(Ae,e)=(D^{-1/2}A^{1/2}A^{1/2}e,D^{1/2}e)=((D^{-1}A)^{1/2}A^{1/2}e,D^{1/2}e)\\ \qquad\qquad \le((D^{-1}A)^{1/2}A^{1/2}e,(D^{-1}A)^{1/2}A^{1/2}e)^{1/2} (D^{1/2}e,D^{1/2}e)^{1/2}\\ \qquad\qquad =(Ae,D^{-1}Ae)^{1/2} (e,De)^{1/2}

ここではシュワルツの不等式を用いた。これを利用して、(Ae,D^{-1}Ae)<<(e,Ae)である場合は(e,Ae)<<(e,De)であることがわかる。

成分を用いてこれを書くと次のとおりになる。

(e,Ae)=\sum_i\sum_ja_{ij}e_ie_j=\sum_i\sum_j\frac{1}{2}\{(-a_{ij})(e_i-e_j)^2+a_{ij}e_i^2+a_{ij}e_j^2\}\\ \qquad\qquad=\sum_i\sum_j\frac{1}{2}(-a_{ij})(e_i-e_j)^2+\sum_i(\sum_ja_{ij})e_i^2\\ \qquad\qquad<<\sum_ia_{ii}e_i^2

ここで、多くの場合に成り立つ\sum_ja_{ij}\simeq0の関係を用いる。また上式は和について成り立つ式であるが、平均的な性質として各iについても成り立つとすると次のとおりになる。

\sum_j\frac{1}{2}|a_{ij}|(e_i-e_j)^2 << a_{ii}e_i^2

つまり以下が成り立つ

\sum_j\frac{|a_{ij}|}{a_{ii}}\frac{(e_i-e_j)^2}{e_i^2} << 1

これはa_{ij}が大きければe_i\simeq e_jが成り立つことを意味している。

つまり、誤差は強い接続上で急激に変化しないという性質を持っていることがわかる。

収束の遅い誤差の例

大きさ1×1の矩形領域における次のような異方性のあるポアソン方程式(Poisson's Equation)

 \{ \begin{array}{ll} \frac{\partial^2\phi}{\partial x^2}+1000\frac{\partial^2\phi}{\partial y^2}=0 \qquad ( \quad in \quad \Omega = \{(x,y):0.0<x,y< 1.0 \} \quad \\ \phi=0 \qquad ( \quad on \quad \partial\Omega \quad ) \end{array}

において乱数を初期値としてGauss Sidel法を用いたSmootherを300反復させた後の誤差の様子である。

 

xの方向よりもyの方向の拡散係数が1000倍大きく、領域の内側の全ての点はyの方向に強い接続、xの方向に弱く接続している。x方向には解がランダムに振動しているのに対してy方向には滑らかに変化している。Smootherで収束の遅い誤差は強い接続上で急激に変化しないということがわかる。


強い接続上で誤差の変化が少ないのでこの強い接続の方向に点を間引き、補間によって間引かれた点の値を表現することができる。このことを利用してF/C SplittiongやProlongationオペレターの作成を行う。

Strong-Connection と Week-Connection

\b{A}\b{e}\simeq 0であった。このとき 行列Aの成分a_{ik}が大きい場合、e_kの値はe_iの値に大きく影響する。逆にa_{ik}が小さい場合、e_kの値がe_iの値に及ぼす影響は小さくなる。

F/C SplitingやProlongationオペレター作成時ではこの強い接続、弱い接続をはっきり区別する必要がある。そこで、i行目の非対角成分a_{ik}\quad(\quad i\ne k\quad)の相対的な大きさから判断して、どの成分kがiに大きく影響するのかを調べ、強い接続、弱い接続に分ける。

\b{A}がM行列(対角が正で非対角が負)の場合は一般的にStrong Connectionは次のようにして定義される。

 -a_{ij}>\theta \max_{k\ne i}\{-a_{ik}\}

\theta(\quad 0<\theta\le1\quad)であるパラメータで0.25とされることが多い。

Prolongation

マルチグリッド法におけるCoarce Grid CorrectionはSmootherが収束しにくいような誤差を落とすためのものであった。

このような収束しにくい誤差の性質として、

ということが分かった。このような誤差を表現できるようなProlongationのオペレータを作ることが目標である

もちろん粗いグリッドの点の値はそのまま反映されるので、i\in Cのときe_i=e_iである。以下ではi\in Fのとき、e_iどのようにして粗いグリッドの点から補間されるのかを調べよう。

Smootherで収束しにくい誤差は残差が小さかった。

(\b{A}\b{e})_i\simeq0

成分ごとに書き出すと次のようになる。

a_{ii}e_i\simeq-\sum_{j\ne i}a_{ij}e_j

つまりe_iの値はiに接続している点全ての値e_jから知ることができる。以下ではこのe_iに接続している点の値について述べる。e_iに接続している点は次の3種類ある。

  1. 粗いグリッドの点j
  2. iに弱く接続している細かいグリッドの点k
  3. iに強く接続している細かいグリッドの点k
細かいグリッドの点iに接続している3種類の点
@:粗いグリッドの点
A:iに弱く接続している細かいグリッドの点
B:iに強く接続している細かいグリッドの点

iに接続している粗いグリッドの点jは値e_jをそのまま使えばいい。

iに弱く接続している細かいグリッドの点k

k\in F^w_iの場合は補間を行わない。次のような値をとることでCoarseGridCorrectionの安定性を増すことができる。
\qquad\qquad\qquad  e_k=e_i

iに強く接続している細かいグリッドの点k

k\in F^s_iの場合は行列の係数の重みを用いた平均として計算する。
\qquad\qquad\qquad  e_k=\frac{\sum_{j\in C_i}a_{kj}e_j}{\sum_{l\in C_i}a_{kl}}
Smootherで収束しにくい成分は強い接続上で変化が少ないのでこのように強い接続間の補完によって値を近似することができる。このような補間はもしもe_jが定数である場合にe_kが同じ値をとるという性質を満たす。このような定数モードの補間ができることは重要である。
分母が0にならないという条件から、全てのiへ強く影響している点kとiと隣接する細かい節点lは必ず一つは粗い節点に隣接していなければならないということがわかる。通常この性質が満たされるようにF/C Splittingを行う

粗いグリッドから細かいグリッドへの補間

節点iの周りの点の値が分かったので節点iの値を求めてみよう。

節点iについての残差が小さいということから

a_{ii}e_i\simeq-\sum_{j\ne i}a_{ij}e_j

が成り立つのであった。周囲の点を上の3パターンに分けると次のようになる。

a_{ii}e_i\simeq-\sum_{j\in C_i}a_{ij}e_j-\sum_{k\in F^s_i}a_{ik}e_k-\sum_{k\in F^w_i}a_{ik}e_k

これに周囲の節点の値を代入して、

a_{ii}e_i\simeq-\sum_{j\in C_i}a_{ij}e_j-\sum_{k\in F^s_i}a_{ik}\frac{\sum_{j\in C_i}a_{kj}e_j}{\sum_{l\in C_i}a_{kl}}-\sum_{k\in F^w_i}a_{ik}e_i

これらを変形してまとめると

(a_{ii}+\sum_{k\in F^w_i}a_{ik})e_i\simeq-\sum_{j\in C_i}\{a_{ij}-\sum_{k\in F^s_i}\(\frac{a_{ik}a_{kj}}{\sum_{l\in C_i}a_{kl}\)\}e_j

e_i\simeq-\sum_{j\in C_i}\(\frac{a_{ij}-\sum_{k\in F^s_i}\frac{a_{ik}a_{kj}}{\sum_{l\in C_i}a_{kl}} }{a_{ii}+\sum_{k\in F^w_i}a_{ik}}\)e_j

AMGにおける粗いグリッドから細かいグリッドへの補間式
w_{ij}=-\frac{a_{ij}-\sum_{k\in F^s_i}\frac{a_{ik}a_{kj}}{\sum_{l\in C_i}a_{kl}} }{a_{ii}+\sum_{k\in F^w_i}a_{ik}}

\b{P}=\[W\\I\]

RistrictionとCoarse Grid Operator

Aが正定値対称行列の場合はRistrictionの行列は次のようにProlongationの転置(Pが複素数行列の場合はHarmetian)を用いるとAノルムの意味での誤差ノルムの最小化ができる。このような方法は有限要素法との類推からGalerkin法とも呼ばれる

R=P^T

A_c=\b{R}\b{A}\b{P}

Aが正定値対称行列ではない場合にはRistrictionにProlongationの転置を用いる方法は必ずしも最適とはいえない。様々なRistrictionを決定する方法が提案されているが、問題依存であることが多い。

F/C Spliting

AMG法ではコースグリッドを行列の係数を元に作成しなければならない。

細かいグリッド上で粗いグリッドの点の集合をC、それ以外の点をFとする。

また細かいグリッドにおいて点iと接続している点j、つまりa_{ik}\ne0の集合をN_iとする。

細かいグリッド上で点iに接続している粗いグリッドの点の集合をC_iというふうに書く。つまりC_i= N_i \cup C

同様に細かいグリッド上で点iに接続している細かいグリッドの点の集合をF_iというふうに書く。つまりF_i= N_i \cup F

粗いグリッドの点から細かいグリッドの点に補間ができるためには次の条件を満たす必要がある。

iを細かいグリッドの点とし、jがiに強く接続している点(つまりj\in S_i)とするとjは粗い点j\in Cか、最低1つは粗い点に強く接続(S_j\cup C\ne \phi)していなければならない

計算効率のことを考えるとできるだけ少ない粗いグリッドの点で細かいグリッドへの補間ができたほうがいいと考えられる。

上の条件を満たしつつできるだけ少ないグリッドで計算を行うためには

粗いグリッド同士は強く接続していない

という性質を満たとよい。この粗いグリッド同士の接続がないという条件はグリッドの数を減らして計算量を小さくするという目的のもので、決定的に重要な条件ではない。

収束について

ここではマルチグリッド法の収束について述べる。解くべき行列Aは正定値対称であるとする。

また、RistrictionI^H_hはProlongationI^h_Hの転置であるとする。つまり、I^H_h=(I^h_H)^T

この場合明らかにコースグリッドのオペレータA^H=I^H_h A I^h_Hも正定値対称である。

ノルム||.||_0は次のようなノルムであるとする。||e||_0=(e,De)、但しDAの対角成分

ノルム||.||_1はオペレーターノルムであるとする。つまり、||e||_1=(e,Ae)

また、ノルム||.||_2を次のように定義する||e||_2=(Ae,D^{-1}Ae)

Prolongationオペレータとコースグリッドコレクションオペレータの値域はオペレータノルムについて直交するという性質について

ProlongationオペレータI^h_HとコースグリッドコレクションオペレータT^hの値域がオペレータノルム||.||_1について直交する、つまりR(I^h_H)\bot_A R(T^h) ということを示す。

v^He^hを任意のベクトルとする。

(I^h_Hv^h,AT^he^h)=(I^h_Hv^h,A\{I-I^h_H(I^H_hAI^h_H)^{-1}I^H_hA\}e^h)=(I^h_Hv^h,Ae^h)-(I^h_Hv^h,AI^h_H(I^H_hAI^h_H)^{-1}I^H_hAe^h)\\ \qquad\qquad\qquad=(I^h_Hv^h,Ae^h)-(v^h,(I^H_hAI^h_H)(I^H_hAI^h_H)^{-1}I^H_hAe^h)=(I^h_Hv^h,Ae^h)-(v^h,I^H_hAe^h)\\ \qquad\qquad\qquad=(I^h_Hv^h,Ae^h)-(I^h_Hv^h,Ae^h)=0


コースグリッドコレクションはプロロンゲーションオペレータの値域にオペレータノルムにおいて直交する部分空間への射影ということがわかる。

以上より

Te=e-I^h_Hv^H

の両辺のAノルムを計算すると

||e||^2_1=||Te||^2_1+||I^h_Hv^H||^2_1

が成り立つ

δを用いた収束の予測

G^mT^mはそれぞれグリッドmにおけるスムーザーとコースグリッドコレクションのオペレータであるとする。
各グリッドmについて以下の式が成り立つとすると
||G^me^m||^2_1\le||e^m||^2_1-\delta||T^me^m||^2_1
ポストスムージングのみのVcycleマルチグリッド法について、誤差のオペレータノルムの収束率が\sqrt{1-\delta}で抑えられる。

以下の式は少し複雑なので、具体的な例を用いて物理的な意味を考えよう。

||G^me^m||^2_1\le||e^m||^2_1-\delta||T^me^m||^2_1

\delta\simeq1つまりマルチグリッド法の収束が早いと仮定すると、

以上からδが1に近いほどコースグリッドコレクションとスムーザーの能力を足し合わせた収束が早いということが言える。

2cycle multigrid 法における証明

上式においてe^hT^he^hを代入すると

||G^hT^he^h||^2_1\le||T^he^h||^2_1-\delta||T^hT^he^h||^2_1

コースグリッドコレクションT^hは射影オペレータであったから、

T^hT^h=T^h||T^he_h||_1\le||e^h||_1が成り立つ。

これを使うと次のように変形できる

||G^hT^he^h||^2_1\le||T^he^h||^2_1-\delta||T^hT^he^h||^2_1=(1-\delta)||T^he^h||^2_1\le(1-\delta)||e^h||^2_1

よって次のようにオペレータノルムでの収束が次のように評価できる。

||G^hT^he^h||_1\le\sqrt{1-\delta}||e^h||_1

δの値について

上のδの式をスムーザーGとコースグリッドコレクションTのオペレータについての2つの分離した式に分けられる。

||G^me^m||^2_1\le||e^m||^2_1-\alpha||e^m||^2_2

||T^me^m||^2_1\le\beta||T^me^m||^2_2
が成り立つ場合は\delta=\frac{\alpha}{\beta}である。

上の式スムーザーGについての式をsmoothing assumptionといい、下のコースグリッドコレクションTについての式をapproximation assumptionという。

2grid Multigrid法での簡単な証明

上の式にe^h=T^h e^hを代入すると

||G^hT^he^h||^2_1\le||T^he^h||^2_1-\alpha||T^he^h||^2_2=||T^he^h||^2_1-\frac{\alpha}{\beta}||T^he^h||^2_1=(1-\frac{\alpha}{\beta})||T^he^h||^2_1\le(1-\frac{\alpha}{\beta})||e^h||^2_1

\deltaを用いた収束の式と見比べると\delta=\frac{\alpha}{\beta}であることがわかる。


\sqrt{1-\delta}がポストスムージングのみのVサイクルマルチグリッド法の収束半径の上限であったので\alphaが大きく、\betaが小さいほど収束が早い

approximation assumptionについて

2cycle multigrid法では、任意のe^hについて \min_{v^H}||e^h-I^h_Hv^H||_0^2\le\beta||e^h||^2_1が成り立つとき
||T^he^h||^2_1\le\beta||T^he^h||^2_2が成り立つ

証明

e^h\in R(T^h)とするとR(I_H^h)\bot_AR(T^h)から、任意のe^Hに対して以下が成り立つ

||e^h||^2_1=(Ae^h,e^h)=(Ae^h,e^h-I^h_He^H)\\ \qquad\qquad = (D^{1/2}D^{-1/2}A^{1/2}A^{1/2}e^h,e^h-I^h_He^H)\\ \qquad\qquad = (D^{-1/2}A^{1/2}A^{1/2}e^h,D^{1/2}(e^h-I^h_He^H))\\ \qquad\qquad = ((D^{-1}A)^{1/2}A^{1/2}e^h,D^{1/2}(e^h-I^h_He^H))\\ \qquad\qquad \le ||(D^{-1}A)^{1/2}A^{1/2}e^h||\quad ||D^{1/2}(e^h-I^h_He^H)||\\ \qquad\qquad = ||e^h||_2\quad ||e^h-I^h_He^H||_0

ここでShwarzの不等式、D^{-1}A>0であることを用いた。

任意の任意のe^hについて \min_{v^H}||e^h-I^h_Hv^H||_0^2\le\beta||e^h||^2_1が成り立つとき、

||e^h||^2_1=||e^h||_2\quad ||e^h-I^h_He^H||_0\le||e^h||_2\sqrt{\beta}||e^h||_1

両辺2乗して

||e^h||^2_1\le\beta||e^h||^2_2

が成り立つ。以上はe^h\in R(T^h)について成り立っていたのでe^h:=T^he^hを代入して最終的に

||T^he^h||^2_1\le\beta||T^he^h||^2_2

が成り立つ。

参考にしたもの

LINK

An Algebraic Mutigrid Tutorial
http://www.llnl.gov/CASC/linear_solvers/talks/AMG_TUT_PPT/sld001.htm

BOOK

Multigrid Methods Stephen F.McCormic 著

PAPER

[1] J.W.Ruge and K.Stuben, algebraic multigrid (AMG), in Mutigrid Methods, S. F. McCormic, ed., Frontiers in Appl. Math. 3, SIAM, Philadelphia, 1987, pp.73-130


Made by Nobuyuki UMETANI  梅谷 信行
n.umetani@gmail.com