数値解析 第12回 (2) ヤコビ法

ヤコビ法のアイデア

 方程式
$\dps{A\xxx=\bbb}$
を成分で書いて、左辺から $a_{ii}\,x_{i}$ の形の項以外を全て右辺へ移してみます: $$ \left\{ \begin{array}{l} a_{11}\,x_1 = b_1 - ( \quad 0\quad\, + a_{12}\,x_2\, + \cdots + a_{1n}\,x_n ) \\ a_{22}\,x_2 = b_2 - ( a_{21}\,x_1\, + \quad 0\quad\, + \cdots + a_{2n}\,x_n ) \\ \qquad\quad\vdots \\ a_{nn}\,x_n = b_n - ( a_{n1}\,x_1 + a_{n2}\,x_2 + \cdots + \quad 0\quad ) \\ \end{array} \right. $$ さらに $a_{ii}$ で割ると $$ \left\{ \begin{array}{l} x_1 = \left\{\, b_1 - ( \quad 0\quad\, + a_{12}\,x_2\, + \cdots + a_{1n}\,x_n ) \,\right\} / a_{11} \\ x_2 = \left\{\, b_2 - ( a_{21}\,x_1\, + \quad 0\quad\, + \cdots + a_{2n}\,x_n ) \,\right\} / a_{22} \\ \qquad\quad\vdots \\ x_n = \left\{\, b_n - ( a_{n1}\,x_1 + a_{n2}\,x_2 + \cdots + \quad 0\quad ) \,\right\} / a_{nn} \\ \end{array} \right. \tag{$\ast$} $$ となります。そこで
アイデア  $(\ast)$ の右辺の $x_{j}$ たちに $\xxx^{(k)}$ の成分を入れたときの、左辺の $x_i$ を $\xxx^{(k+1)}$ の成分としよう。

実行例

 試しにやってみましょう。
Ex.1  $ \dps{ \left\{ \begin{array}{l} 3x-2y = 1\\ \phantom{3}x+3y = 4\\ \end{array} \right. } $
 $(\ast)$ は
$\dps{ \left\{ \begin{array}{l} \mbox{新} x = \dps{\frac{1}{3}} \left(1+2\times\mbox{旧}y\right) \\ \mbox{新} y = \dps{\frac{1}{3}} \left(4-\mbox{旧}x\right) \\ \end{array} \right. }$
となり、初期値 $(x,y)=(0,0)$ で実行すると
Step 0 : (x,y) = (0.0000000000, 0.0000000000)
Step 1 : (x,y) = (0.3333333333, 1.3333333333)
Step 2 : (x,y) = (1.2222222222, 1.2222222222)
Step 3 : (x,y) = (1.1481481481, 0.9259259259)
Step 4 : (x,y) = (0.9506172840, 0.9506172840)
Step 5 : (x,y) = (0.9670781893, 1.0164609053)
Step 6 : (x,y) = (1.0109739369, 1.0109739369)
Step 7 : (x,y) = (1.0073159579, 0.9963420210)
Step 8 : (x,y) = (0.9975613474, 0.9975613474)
Step 9 : (x,y) = (0.9983742316, 1.0008128842)
Step10 : (x,y) = (1.0005419228, 1.0005419228)
Step11 : (x,y) = (1.0003612819, 0.9998193591)
Step12 : (x,y) = (0.9998795727, 0.9998795727)
Step13 : (x,y) = (0.9999197151, 1.0000401424)
Step14 : (x,y) = (1.0000267616, 1.0000267616)
Step15 : (x,y) = (1.0000178411, 0.9999910795)
Step16 : (x,y) = (0.9999940530, 0.9999940530)
Step17 : (x,y) = (0.9999960353, 1.0000019823)
Step18 : (x,y) = (1.0000013216, 1.0000013216)
Step19 : (x,y) = (1.0000008810, 0.9999995595)
Step20 : (x,y) = (0.9999997063, 0.9999997063)
Step21 : (x,y) = (0.9999998042, 1.0000000979)
Step22 : (x,y) = (1.0000000653, 1.0000000653)
Step23 : (x,y) = (1.0000000435, 0.9999999782)
Step24 : (x,y) = (0.9999999855, 0.9999999855)
Step25 : (x,y) = (0.9999999903, 1.0000000048)
Step26 : (x,y) = (1.0000000032, 1.0000000032)
Step27 : (x,y) = (1.0000000021, 0.9999999989)
Step28 : (x,y) = (0.9999999993, 0.9999999993)
Step29 : (x,y) = (0.9999999995, 1.0000000002)
Step30 : (x,y) = (1.0000000002, 1.0000000002)
Step31 : (x,y) = (1.0000000001, 0.9999999999)
Step32 : (x,y) = (1.0000000000, 1.0000000000)
めでたく厳密解 $(x,y)=(1,1)$ に収束しました。

ヤコビ法のアルゴリズム

$\xxx^{(k)}$ の成分を $$ \xxx^{(k)}=\mat{c}{x^{(k)}_1 \\ \vdots \\ x^{(k)}_n}$$ と書けば、アルゴリズムは次の通りです。
Algorithm 2 ( ヤコビ法 ) 
  1. $\xxx^{(0)} = $ 適当なベクトル, $k=0$
  2. for $i=1$ to $n$
       $x^{(k+1)}_i = \left.\left\{ b_i-\dps{\sum_{j=1}^{i-1}a_{ij}\,x^{(k)}_j} -\dps{\sum_{j=i+1}^na_{ij}\,x^{(k)}_j} \right\}\right/a_{ii}$,
    $k=k+1$
  3. if ( 終了条件 ) $\xxx^{(k)}$ を出力して終了 else 2° へ戻れ
 終了条件としては次のような条件を設定します。
  • $\xxx^{(k+1)}-\xxx^{(k)}$ の成分が全て十分に小さい、または
  • 反復回数が十分に大きい
反復法の常で、収束が保証されていないのでループ回数に上限を設定します。