2004/07/06 先週はヤコビ法とガウス=ザイデル法を学んだ # その原理は.. ? ## まず、方程式系を一般化して考える f_1(x_1,x_2,..,b_1) = 0 f_2(x_1,x_2,..,b_2) = 0 .. f_n(x_1,x_2,..,b_n) = 0 を、 x_1 = g_1( x2,x3,..,x_n,b_1) x_2 = g_2(x1, ,x3,..,x_n,b_1) .. x_n = g_n(x1, ,x3,.., b_1) と変形する。 g_i を素直に使うのがヤコビ法 一つ前の x_i を、次の x_{i+1} の計算に利用するのが、ガウス=ザイデル法 値が収束する場合は、収束点の廻りをいったり来りする ( 計算の無駄 ?? ) 値の収束を、片側からだんだんと収束させる ( 数値的に安定する ?? ) SOR 法 ( Successive Over-Relaxation Method : 逐次過剰緩和法 ) SOR 法の原理 ガウス=ザイデル法での修正をα倍して修正する 0 < α < 2 ( α == 1 の時がガウスザイデル法と一致する ) # αの範囲が上記の範囲であることが、SOR 法が収束するた # めの必要条件であることがわかっている。 ## 証明は、オストロフスキーの定理 ### 仮定 : A が対称で、正定値の場合 ## この定理は、どのαが高速かを示してくれない 経験的には、α = 1.8 から 1.99 を利用するとよいことが判っている SOR 法 \delta x_i : 修正量 この修正量は、例えば、ガウス=ザイデル法で求める {}_{k+1}x_i = {}_kx_i + \alpha \delta x_i # 具体的には、 r_i = b_i - \Sum_{j=1}^n a_{ij} x_j {}_{k+1}x_i = {}_kx_i + \alpha + \alpha \frac{r_i}{a_{ij}} # 計算式は、テキストのガウス=ザイデル法と同じ ## 前回、テキストの方法は、安定度が落る ( 誤差を含み易 ## い.. ) と述べたが、今回の場合は、問題が生じない。 ### これは、前回が値そのものの計算だったのに、今回は残差の計算だから {}_{k+1}x_i = \frac{b_i - \Sum_{j=1}^{k-1} a_{ij}x_j - \Sum_{j=k+1}^{n} a_{ij}x_j}{a_{ii}} \delta_k x_i = {}_{k+1}x_i - {}_{k}x_i {}_{k+1}\~x_i = {}_k x_i + \alpha \delta_k x_i # SOR 法は、流体の計算に利用される ## ヤコビ、ガウス=ザイデル、SOR 法では SOR 法が一番よく利用されている ### 収束も速いし、経験的にも判っているので..、性質もよい ### cf. バーガ著、渋谷政昭訳、「大型行列の反復解法」、サイエンス社, 1972 ### 30 年前の本でも役に立つ # 飛行の羽の解析 # 羽の表面が問題 # 都合がわるければい、αを 1 にして、ガウス=ザイデルにすればよい。 == 解法の分類 直接法 ガウスの消去法できまり 間接法 古典的 SOR 法が良く利用されている 最近 CG 法 # これらの方法は、100 年前 ( つまり、計算機がなかったころ.. ) から使われている # 収束の方法が目に見える # 一度でも計算回数を減らしたい # 計算回数を減らすための工夫が色々されている CG 法 ( Xonjugate Gradient 法 : 共役勾配法 ) 元の方程式に予想値を代入すると、残差が出る 残差を 0 にするということが、方程式を解くことになる 残差に関する評価関数 ( 恒に非負値を取り、残差が 0 の時に 0 [最小値] となる ) を 考え、評価関数の最小値を求めるということを考える => 最適値計算の手法が利用できるようになる [例] 今いる場所から、できるだけ、山の頂点に行きたい場合 => できるだけ上りの傾斜がきつい方向に行けばよい # 山上り法 n 次元の場合は、直交性のある操作で n 回で解に到達する方法があるとよい => CG 法 # 丸め誤差がなければ、n 回で解に到達することが証明されている ## しかし、実際には、丸め誤差が必ず入るので、解けない問題が出てくる ### 「直交性」を利用しているがこの「直交性」が丸め誤差で崩れてしまう # (n 回なので..) 一時はもてはやされたが、(誤差があるので) すたれた CG がうまくゆく場合ゆかない場合 評価関数の断面の形が円に近い.. うまくゆく 各方向に同じような性質がある 固有値の大きさが同じような大きさ 評価関数の断面の形がつぶれている 方向によって、異る性質をもつ 固有値の大きさが桁違いに異る 問題 ( 行列 ) を変形 ( 前処理 ) して、CG 法で解くのに都合の良い形に変形して 解く 前処理 CG 法 # 与えられた行列 A を、固有値が同じ大きさになる行列 B に変形 ## 変形は問題 ( 行列 A ) 毎に異ることに注意 ## 変形そのものが難しい問題になっていることも多い CG 法の振舞い 各次元毎に、残差を 0 にする。 一つの次元を 0 にする操作が、残差のノルムを減らすという保証がない 計算の途中で、一時的に残差ノルムが増大する可能性もある # ほとんどの場合は、残差ノルムも減るのだが.. ## n 回目には必ず残差ノルムが 0 になることは保証されている CG 法でかつ、残差が減ることを保証した方法がある ファン・デ・ア・フォルスト(オランダ) の方法 # 現在は、前処理をしたファン・デ・ア・フォルスト法の CG 法が一番安心して使える ## この方法は、まだ発展段階.. CG 法の前処理をどうするか ? A -> B の変形 前処理の極端な例 A * I -> B = A 何もしない A * A^{-1} -> B = I すでに解けている # 「前処理に時間をかければ後処理で楽ができる」という原理 # 一般には、I と A^{-1} の間の行列をかけることになる。 == 応用 流体 ガウスの消去法 SOR 法 構造計算 CG 法 ( 前処理 + ファン・デ・ア・フォルスト法 ) == 最近の手法 (CG 法) クリロス空間法 まだ、良いかどうかわからない # 新しい方法がでるということは、以前の方法が不十分ということ == 今週の課題 SOR 法 r_i = b_i - \Sum_{j=1}^{n} a_{ij} x_i {}_{k+1}x_i = {}_{k}x_i + \alpha \frac{r_i}{a_{ii}} \alpha = 1.80 から 1.99 1.9 で試す 余裕がある人は、他の値も試す 入力データは、先週のガウス=ザイデル法と同じ物 == 前期の講義は今日で御仕舞い 来週 まとめ ( 9:00 - 9:45 ) テスト ( 9:45 - 10:45 ) # テストが終った人が、外で騒ぐので、他の講義に迷惑がかかる o 山は掛けよう 必ず当る 「問題」を、去年の人に尋ねよう 「解答」は、去年の人にきいてはだめ ( 答が当っているとい保証はない ) # 落とすテストではない、将来に役立つような知識を得てもらうことが目的 # 「まとめ」の内容は、答をしゃべっているようなもの