SOR法

SOR法︵SORほう、英: Successive Over-Relaxation、逐次加速緩和法︶とは $n$$n$元連立一次方程式$A{\boldsymbol {x}}={\boldsymbol {b}}$$A{\boldsymbol {x}}={\boldsymbol {b}}$を 反復法で解く手法の一つであり、 ガウス＝ザイデル法に加速パラメータ$\omega$$\omega$を導入してその修正量を拡大することで、 更なる加速を図った手法である。

反復のスキーム

$n$$n$次正方行列$A$$A$は、上三角行列$U$$U$、 下三角行列$L$$L$、対角行列$D$$D$の和に分離すると、

$A=L+D+U$$A=L+D+U$
と書ける。

非対角成分に相当する項をすべて右辺に移項し、すべての量$x_{1},x_{2},\ldots ,x_{n}$$x_{1},x_{2},\ldots ,x_{n}$に 各段階で得られている最新のデータを代入するようにする(ガウス＝ザイデル法)。こうして計算された値を${\boldsymbol {\tilde {x}}}^{(k+1)}$${\boldsymbol {\tilde {x}}}^{(k+1)}$ とすると、${\tilde {x}}_{i}^{(k+1)}$${\tilde {x}}_{i}^{(k+1)}$は次の形となる。

${\displaystyle {\tilde {x}}_{i}^{(k+1)}={\frac {1}{a_{ii}}}\left(b_{i}-\sum _{j=1}^{i-1}a_{ij}x_{j}^{(k+1)}-\sum _{j=i+1}^{n}a_{ij}x_{j}^{($k)}\right)} ${\tilde {x}}_{i}^{(k+1)}={\frac {1}{a_{ii}}}\left(b_{i}-\sum _{j=1}^{i-1}a_{ij}x_{j}^{(k+1)}-\sum _{j=i+1}^{n}a_{ij}x_{j}^{(k)}\right)$
この値を次段でそのまま採用せずに、ガウス＝ザイデル法で本来修正される量${\displaystyle {\tilde {x}}_{i}^{(k+1)}-x_{i}^{($k)}} ${\tilde {x}}_{i}^{(k+1)}-x_{i}^{(k)}$に1より大きい 加速パラメータ︵では緩和因子、緩和係数と呼ばれている︶$\omega$$\omega$を乗じてこの修正量を拡大し、これを前段の近似値${\displaystyle x_{i}^{($k)}} $x_{i}^{(k)}$に加えることで、新たな値は

${\displaystyle x_{i}^{(k+1)}=x_{i}^{($k)}+\omega \left({\tilde {x}}_{i}^{(k+1)}-x_{i}^{(k)}\right)} $x_{i}^{(k+1)}=x_{i}^{(k)}+\omega \left({\tilde {x}}_{i}^{(k+1)}-x_{i}^{(k)}\right)$
とできる。ただし、桁落ちを防ぐ観点からこの式の通り計算するのではなく、

${\displaystyle x_{i}^{(k+1)}=(1-\omega )x_{i}^{($k)}+\omega {\tilde {x}}_{i}^{(k+1)}} $x_{i}^{(k+1)}=(1-\omega )x_{i}^{(k)}+\omega {\tilde {x}}_{i}^{(k+1)}$
として計算するか、または本節の最後に書かれた式を用いるのがよい。

この漸化式を、上の$A=L+D+U$$A=L+D+U$を用いて行列で表現すると、

${\displaystyle {\Biggl \{}{\begin{matrix}{\boldsymbol {\tilde {x}}}^{(k+1)}&=&D^{-1}\left({\boldsymbol {b}}-L{\boldsymbol {x}}^{(k+1)}-U{\boldsymbol {x}}^{($k)}\right)\\{\boldsymbol {x}}^{(k+1)}&=&{\boldsymbol {x}}^{(k)}+\omega \left({\boldsymbol {\tilde {x}}}^{(k+1)}-{\boldsymbol {x}}^{(k)}\right)\end{matrix}}} ${\Biggl \{}{\begin{matrix}{\boldsymbol {\tilde {x}}}^{(k+1)}&=&D^{-1}\left({\boldsymbol {b}}-L{\boldsymbol {x}}^{(k+1)}-U{\boldsymbol {x}}^{(k)}\right)\\{\boldsymbol {x}}^{(k+1)}&=&{\boldsymbol {x}}^{(k)}+\omega \left({\boldsymbol {\tilde {x}}}^{(k+1)}-{\boldsymbol {x}}^{(k)}\right)\end{matrix}}$
となり、この2式から${\boldsymbol {\tilde {x}}}^{(k+1)}$${\boldsymbol {\tilde {x}}}^{(k+1)}$を消去することで、次式が得られる。

${\displaystyle {\boldsymbol {x}}^{(k+1)}=\left(I+\omega D^{-1}L\right)^{-1}\left\{\left(1-\omega \right)I-\omega D^{-1}U\right\}{\boldsymbol {x}}^{($k)}+\omega \left(D+\omega L\right)^{-1}{\boldsymbol {b}}} ${\boldsymbol {x}}^{(k+1)}=\left(I+\omega D^{-1}L\right)^{-1}\left\{\left(1-\omega \right)I-\omega D^{-1}U\right\}{\boldsymbol {x}}^{(k)}+\omega \left(D+\omega L\right)^{-1}{\boldsymbol {b}}$
上式における${\displaystyle {\boldsymbol {x}}^{($k)}} ${\boldsymbol {x}}^{(k)}$の係数 $M_{\omega }=\left(I+\omega D^{-1}L\right)^{-1}\left\{\left(1-\omega \right)I-\omega D^{-1}U\right\}$$M_{\omega }=\left(I+\omega D^{-1}L\right)^{-1}\left\{\left(1-\omega \right)I-\omega D^{-1}U\right\}$を反復行列という。

実際の数値計算においては、これを各成分について表した下の式が用いられる。

$x_{i}^{(k+1)}=(1-\omega )x_{i}^{(k)}+\omega {\frac {1}{a_{ii}}}\left(b_{i}-\sum _{j=1}^{i-1}a_{ij}x_{j}^{(k+1)}-\sum _{j=i}^{n}a_{ij}x_{j}^{(k)}\right)$

収束性

反復行列の固有値を$\lambda$$\lambda$とすると、

$\max |\lambda |\geq |\omega -1|\quad \forall \omega$$\max |\lambda |\geq |\omega -1|\quad \forall \omega$
が成立することから、少なくとも$0<\omega <2$$0<\omega <2$でなければSOR法の収束性は保証されない 。

さらに、正定値対称行列$A$$A$を係数にもつ方程式$A{\boldsymbol {x}}={\boldsymbol {b}}$$A{\boldsymbol {x}}={\boldsymbol {b}}$に対するSOR法は、 加速パラメータ$\omega$$\omega$が$0<\omega <2$$0<\omega <2$のとき必ず収束する︵Ostrowskiの定理︶。

また、$\omega =1$$\omega =1$のときガウス＝ザイデル法と同じになり、$\omega$$\omega$が1より小さいとき ガウス＝ザイデル法より収束が遅くなる。ただし、ガウス＝ザイデル法で収束しないような問題には使える 。

加速パラメータ$\omega$の選択

一般に加速パラメータ$\omega$$\omega$の値をあらかじめ最適に定めることはできない。そのため、問題ごとに適当な値を選択する必要がある。

しかし、$\omega$$\omega$の最適な値を決定することができる例も存在する。それは、係数行列$A$$A$が、ある基本行列$P$$P$に 対して

$PAP^{-1}=\left[{\begin{matrix}D_{1}&M_{1}\\M_{2}&D_{2}\end{matrix}}\right]$$PAP^{-1}=\left[{\begin{matrix}D_{1}&M_{1}\\M_{2}&D_{2}\end{matrix}}\right]$
という形の行列に相似変換することができ、さらにヤコビ法の反復行列$M_{J}=-D^{-1}\left(L+U\right)$$M_{J}=-D^{-1}\left(L+U\right)$の スペクトル半径$\rho \left(M_{J}\right)$$\rho \left(M_{J}\right)$が既知であるときである。 なお、上の行列内の$D_{1},D_{2}$$D_{1},D_{2}$は対角行列である。

このとき、SOR法の反復行列のスペクトル半径$\rho \left(M_{\omega }\right)$$\rho \left(M_{\omega }\right)$が最小となる$\omega$$\omega$の最適値は、 次の形で得られる。

$\omega ={\frac {2}{1+{\sqrt {1-\rho \left(M_{J}\right)^{2}}}}}$

近年の研究

共役勾配法をはじめとしたクリロフ部分空間法の普及が進んだことでSOR法の使用が減ってしまったこともあったが、離散勾配法 (構造保存型数値解法の一つ) との関係が明らかになったことで再び注目されている。

脚注

(一)^ ^abcdefg山本哲朗﹃数値解析入門﹄︵増訂版︶サイエンス社︿サイエンスライブラリ 現代数学への入門14﹀、2003年6月。ISBN 4-7819-1038-6。 

(二)^ 幸谷智紀 (13 October 2007). 連立一次方程式の解法2 -- 反復法 (PDF) (Report). 2018年3月30日閲覧。

(三)^ “SOR法” (2004年12月16日). 2018年3月30日閲覧。

(四)^ Varga, R. S. (2009). Matrix iterative analysis (Vol. 27). Springer Science & Business Media.

(五)^ 宮武勇登, 曽我部知広, & 張紹良. (2017). 微分方程式に対する離散勾配法に基づく線形方程式の数値解法の生成. 日本応用数理学会論文誌, 27(3), 239-249.

(六)^ Miyatake, Y., Sogabe, T., & Zhang, S. L. (2018). On the equivalence between SOR-type methods for linear systems and the discrete gradient methods for gradient systems. en:Journal of Computational and Applied Mathematics, 342, 58-69.

参考文献

●森正武﹃数値解析﹄共立出版、2002年2月。ISBN 4-320-01701-3。 

●杉原正顯, 室田一雄, 線形計算の数理, 岩波書店, 2009年.

●線形方程式の反復解法、 一般社団法人 日本計算工学会 編、藤野清次 著、阿部邦美 著、杉原正顯 著、丸善出版、2013年09月。

●数値線形代数の数理とHPC, 櫻井鉄也, 松尾宇泰, 片桐孝洋編︵シリーズ応用数理 / 日本応用数理学会監修, 第6巻︶共立出版, 2018.8

関連項目

• 数値解析
• 数値線形代数
• ヤコビ法
• ガウス＝ザイデル法