2004/06/29 直接法で連立方程式を解く場合は、Gauss 法が最適 ただし、条件がある。 一般に、方程式を解く方法は、二種類ある 直接法 ( Gauss 法 etc.. ) 反復法 解の推定が可能ならば、それを方程式に放り込み修正を行う方法 方程式を解く場合に一般的に適用できる方法 直接法があるのも拘らず、反復法を利用する場合がある理由は ? 問題が小さい ( i.e. 演習で解く問題など.. ) 場合は、直接法で十分 直接法は、密行列を対象としている。 現実の問題では.. ? 大規模偏微分方程式を解く 有限差分法 .. スパース行列 有限要素法 .. スパース行列 境界要素法 .. 密行列 有限要素法 1 次元 ( 3 重対角 ) \\ 0 \\\ \\\ \\\ 0 \\ 2 次元 ( 5 重対角 ) \\\ \ 0 \\\ 0 \ \ \\\ \ \ 0 \\\ 0 \ \\\ 3 次元 ( 7 重対角 ) \\\ \ \ 0 \\\ \ 0\ \ 0 \\\ 0 \ \ 0\ \\\ 0 \ \ \\\ 一般には、有限差分法、有限要素法、つまりスパース行列を解くことが多い スパース行列では、0 の要素が多い 直接法を使うと、せっかく元々 0 だった要素が計算の途中で非0 になる ( fill-in という現象 ) 0 のままなら計算はサボれるがそうはいかない スパースでもあまり嬉しくない 反復法では、行列の要素が 0 だと、計算しなくてすむ # 本質的に Ax を計算することなので.. 密行列 スパース行列 消去法 ○ × 反復法 × ○ 複数の問題を解く場合 直接法 問題が変っても計算量は変らない 反復法 推定される答の候補が解に近ければ速く計算できる 発展的な問題 時刻によって、解が変ってゆく問題 時刻の間隔が短ければ解も近い # 物理現象は保存則が成立するので.. 前の時刻の解を次の時刻の解の初期値の候補として使えば収束が速くなる可能性が高い 反復法が有利な問題 大きな次元のスパース行列を扱う場合 時間発展問題を解く場合 直接法の利点 ( 反復法の欠点 ) 初期値が不要 # ここまでは前置き == ガウス＝ザイデル法 与えられた問題 Ax = b を書き下すと.. a_11 x_1 + a_12 x_2 + a_13 x_3 = b_1 a_21 x_1 + a_22 x_2 + a_23 x_3 = b_2 a_31 x_1 + a_32 x_2 + a_33 x_3 = b_3 # 対角要素に着目して書き直すと.. x_1 = \frac{1}{a_11}(b_1 - a_12 x_2 - a_13 x_3) x_2 = \frac{1}{a_21}(b_2 - a_21 x_1 - a_23 x_3) x_3 = \frac{1}{a_31}(b_3 - a_31 x_1 - a_22 x_2 ) # これを f(v) ととして.. # ここで、A は密行列として式を立てているが、スパースの場合は、 # 要素が 0 の部分は計算を省略することができるので計算が略 {}_0x_1 {}_0v = ( {}_0x_2 ) {}_0x_3 として、v_{i+1} = f(v_{i}) と計算する (ヤコビ法)。 {}_1x_1 = \frac{1}{a_11}(b_1 - a_12 {}_0x_2 - a_13 {}_0x_3) {}_1x_2 = \frac{1}{a_21}(b_2 - a_21 {}_0x_1 - a_23 {}_0x_3) {}_1x_3 = \frac{1}{a_31}(b_3 - a_31 {}_0x_1 - a_22 {}_0x_2 ) # これは、元の行列の性質がよければ収束し、解になる。 # 収束の精度はそれほど気にしなくてもよい。 ## 国土地理院では、精度が欲しい ## 応用では長い距離の差分が必要なのだが、差分を取ると精度が落るので予め精度が十分に欲しい # ヤコビ法を初心者にプログラミングさせると、多くの場合自然に、 # ガウス＝ザイデル法になってしまう。 ## Why ? 収束する変数を新旧で共用してしまうと自動的にガウスザイデル法になる cf. Text p.55 ガウスザイデル法 ( 新旧の変数を同じメモリで共有する ) {}_1x_1 = \frac{1}{a_11}(b_1 - a_12 {}_0x_2 - a_13 {}_0x_3) {}_1x_2 = \frac{1}{a_21}(b_2 - a_21 {}_1x_1 - a_23 {}_0x_3) {}_1x_3 = \frac{1}{a_31}(b_3 - a_31 {}_1x_1 - a_22 {}_1x_2 ) # 実は、ガウス＝ザイデル法の方がヤコビ法より収束性が高いので一般 # の応用ではこちらを使う事が多い ## 殆どの場合、ヤコビ法よりガウス＝ザイデル法の方が良いが並列化を行う場合は、ヤコビ法の方がよい。 対角優位ならば収束する。 普通の偏微分方程式なら、 他の要素は 1 対角要素は 4 + h が成立するので、普通、この性質が成立する # 極端な話、対角の値が 0 だと、zero div で止る ## 対角優位は、十分条件。一般には、対角要素の値が、他の要素より大きければ、大丈夫なことが多い。 物が壊れない時は、それを表す行列は、安定 (対角優位) な行列になっている 不安定な現象の場合 ( cf. プラズマ現象 ) は、不安定な行列がでることもある # 不安定な例 ( 電気回路 ) # # +---------------------------------- # | # _ # - (電池) # | # +---------------------------------- # (無限に長い抵抗) # 安定度はぎりぎり ガウス＝ザイデル法を使う場合は、安定しているかどうかが問題 物理現象が安定しているならば、使ってよい 不安定な物は.. ?? # 物を作る場合は、安定なものしか使っては行けない。 # ハンドドリルで穴を空ける場合 # まるい穴はできない、オニギリ型になる # ドリルの三つの歯の一つがひっかかって他の歯が梃の仕組でけずる # どの穴も異る形の穴になる ( 不安定 ) # == 予定 6/29 本日 ( ガウス=ザイデル法 ) 7/6 講義 SOR 法 7/13 まとめ + 試験 ( 1 時間位 ) == 本日の演習 [基本] ガウス＝ザイデル法のプログラムを作成し提出する [応用] ヤコビ法のプログラムも作り、両者で同じ問題を解き、収束の差を調べる # 二つの結果を比較しよう !! Text の図 3.12 は、ガウス＝ザイデル法 Text の Program の説明 # 内側の Loop S = b[i] for ( j = 0; j < n; j++ ) { S -= a[i][j]*x[j] S = S/a[i][i] } x[i] += S; # S に a_11 x[i] を加えて、また、a_11 x[i] を引いている ## 数学的には正しい ## 数値的には、問題かも !! ### a_11 x[i] の絶対値が大きい ( 対角優位だと大きくなる !! ) と誤差が増える。 # 嫌な人は、Loop を分けよう !! 問題の係数は Text の p. 58 の物を使う。答は p.57 にある == 最近の流行 共役傾斜法 CG 法 # 一般の反復法は、繰り返し回数を予測することはできない # CG 法は、無限桁で計算すれば、n 回で収束するという証明がある # 現実には、無限桁では計算できないので、一旦廃れた # 最近、この問題を回避することによって、n 回に近い回数で収束するようになり見直された ## 直接法は、ガウス法で、決りだが、反復法は、まだまだ「一番」という答がない ## 研究するならば、反復法 !!