2005/09/27 数値解析学と演習 福井 先生 text p.267 今回と次回は、固有値問題のヤコビ法を学ぶ # 前期後期の内容で、ヤコビ法が最もややっこしい。 # 線型方程式は難しくない、後期の積分法などは、補間法応用なので理解しやすい [定義] 固有値 与えられた行列 A に対して A x = \lambda x を満す \lambda の事 ( cf. 線型代数の教科書 ) # x は、\lambda に対する「固有ベクトル」と呼ぶ # 何が要求されているのだろうか ? 右辺での A は、 n^2 の情報がある、しかし、左辺は n の情報しかない.. => 一種の制約条件を求めている。 [例] 紐は、そのままだと、その両端は、かなり自由に動ける 一方の端点を固定しただけでも、比較的自由 => 両端点を固定すると、かなり不自由 # 両端点を壁に止めるのは大変だが.. 両端点を固定したときの運動のようす 両端点を固定されているので、その紐の長さの 1/n の周波数を持つ振動しかできない => 固有振動数 ( 特定の限られたモードを調べる ) 建物の固有振動数を求めることが、固有値問題のよくある応用例 # 数学的には、 Ax = \lambda x の形を取る問題は多いが.. # 応用としては、この形が多い 最近の地震で、高いビル ( cf. 六本木ヒルズ.. ) のエレベーターが壊れた 低いビルのエレベーターは大丈夫だったのだが.. => 地震の波と建物の固有振動数に近かったので # 建物が地震に影響されるのは、地震の振動数と建物の固有振動数が近い場合 !! # 共振の応用の例 ラジオのチューニング ( バリコンを変えると共振して繋がる ) ( cf. バリコン : バリアブルコンデンサー ) == [どうやって解くか ?] cf. 線型方程式の場合 A x = b ---(変形)---> A' x = b' a_11 a_12 a_13 A = ( a_21 a_22 a_23 ) a_31 a_32 a_33 a_11' 0 0 A' = ( a_22' 0 ) 0 0 a_33' # 同様にして.. A x = \labda x ---(変形)---> A' x = \lambda x # 要するに、「対角化」すれば良い # 問題は.. 「答を変えない行列の操作」が「何か」という問題 # => 線型方程式では、「基本操作」がそれ 固有値を変えない行列の変形操作は何か ? # 固有値を変えてしまう操作では違う問題を解いてしまうことになる !! # cf. 基本変形だと、固有値が変る => 代数方程式とは異る操作が必要 !! 「相似変換」 ( cf. text p.269 ) 固有値を変えない行列の変形操作 相似変換とは ? 大きさ (長さ ) を変えない、座標変換 # ゲームの中で、「絵がクルクル回る」のは、(回転なので..) 相似変換を行っている例 # i386 の SSE (座標変換機能) も、固有値問題を解くのに使える # # 精度(桁数)が足りないのが問題だが.. 2 x 2 の固有変換 2 個の固有値 # 一般に n x n の相似値は、n 個ある 単位円を潰して回転させた楕円 楕円の軸 ( 長軸と短軸の二つある ) の長さが固有値 相似変換 楕円の大きさと形を変えずに、原点を中心にして回転し、長軸と、短軸が、x, y 軸に重なるようにする => すると、軸の長さ ( = 固有値 ) が簡単に計算できる。 # [注意] 通常、固有値が異るが、同じ場合は、円になる # 三次元の場合は、ラグビーボールの潰れたもの ## 四次元の場合 ... => 各自、想像してみよう.. == [固有値問題と線型方程式との違い] 線型方程式では、一旦 0 にした要素は、0 のまま # もし、手順を間違えなければだが... # => 必ず、一定の回数で終了することを示すことができる <=> 固有値問題の場合は、一旦 0 にした要素が、0 でなくなる可能性がある => 繰り返し回数が一定でなくなってしまう !! # これは、どんな意味をもつのか ? 固有値を解く ( 数学的な.. ) 方法 A x = \lambda x より、 (A - \lambda) x = 0 従って、 |A - \lambda| = 0 ここで、左辺は、 \lambda に関する n 次多項式になっている。 固有値問題を解く事と n 次元の代数方程式を解くのは同じ # 現在は、色々な(数値的に..)問題があるので取られないが... 昔 ( まだ、現在の良い方法が知られていなかった時代.. ) は、固有値問題を、代数方程式を解くことで解いていた。 n が 4 以下ならば、「公式」が存在する => つまり、「一定の回数で解くことができる」ということ n が 5 以上 => 公式が存在しない ( アーベルの定理 ) => 一定回数で解法が存在しない # 代数方程式の解法は、前期に学んだ == [固有値の数値解法] ヤコビ法 ( 対称行列の固有値問題の解法 ) 物理現象を ( 普通に.. ) 定式化すれば、対称行列になる # 対称性は、「物理の『保存則』」に相当する # 固有値問題では、固有値が異る方が良い => 同じだと計算が面倒 固有値が異れば、固有ベクトルは互いに直交している # 二次元では、円になる 固有値が同じだと固有ベクトルが決らない # 実際には、固有ベクトルに自由度が現れるので扱いにくい # もし、問題そのものに対称性があると、この結果 ( 保存則が原因でない.. ) 対称性が生まれ、固有値が重複してしまう # => これは計算の上で都合が悪いので、問題を解く ( 定式化する ) 前に、その問題自身の対称性は取り除く !! # 「正方形の面」 # => 対称軸が沢山 ( 4 つある ) # 対称性質を取り除いた 1/8 の三角形で考えればよい !! [相似変換の数学的な定式化] M : 正則とする M (A x) = M (\lambda x) = \lambda M x より、 y = M x と置けば x = M^(-1) y なので、 M A M^(-1) y = \lambda y つまり、 A => M A M^(-1) は、固有値を変えなり変換となっている ( 相似変換 ) # 一般に、M が正則であれば、どのような場合でも相似変換となるが、 # もっと簡単な M に関してのみ扱う ヤコビ法で利用する M_pq(\thita) M_pq(\thita) = { m_ij } m_ji = \cos{\thita} i = j = p | i = j = q \sin{\thita} i = p, j = q - \sin{\thita} i = q, j = p 0 i \ne j 1 i = j これは、p 次元目と q 次元目が作る平面上の \thita 回転になっている M の性質 M^{-1} = M^t # 逆行列は、転置行列になっている # 逆行列を計算する必要がないので大変嬉しい # # というか、必要がないような M を昔の人が考えた # # # 逆行列の計算が必要だったら、(計算機もないわけだから.. ) そもそもできなかっただろう.. M_pq を使えば、A の p,q 要素を 0 にできる しかし、二つの行と二つの列に影響を及ぼすので、一旦 0 にした部分が、また、0 でなくなってしまう可能性がある !! => 三重対角行列までは、問題なく、0 にできるが、対角要素の前後左右の要素まで 0 にしようとすると、他の場所にも影響を及ぼしてしまい、他の要素に 0 でない要素ができてしまう。 [数値解法] 数学的には、三重対角までしか変形できない => つまり「解けない」ってこと ??? Yes : 実際に公式もないわけだし... # 逆に解けるとすれば、「公式がある」ことを意味する。 => じゃあ.. ? ところが、ヤコビ法は一つ「良い」性質がある 対角要素以外は、絶対値が小くなる # 他の要素のエネルギーが、対角要素に移動する !! => 数学的には 0 にできないが、数値的には、「0 に近く」できるので、対角要素以外の要素の大きさが十分に 0 に近くなったら 0 とみなしても問題ない !! ヤコビ法で、計算を速くするには、 絶対値が最も大きい物を 0 にすればよい !! [具体的な方法] text p.270 の (7) 式を利用して、p,q 成分を 0 にするような \thita を求める \thita (根は二つあることに注意) は次の方程式を満すことになる \cos{\thita}\sin{\thita}(a_qq-a_pp)+(\cos^2{\thita}-\sin^2{\thita}a_pq = 0 ここで、 \thita = \frac{1}{2}\arctan{\frac{2 a_pq}{a_pp-qa_qq} となる。 (prof) 今、a_pq \ne 0 とする ( a_pq = 0 の場合は \thita を求める必要がない ) t = \tan{\thita} とすると、元の方程式の両辺を \cos^2{\thita} で割れば、 t ( a_pp - a_qq ) + a_pq ( 1 - t^2 ) = 0 となるので、これより、二次方程式を解けば t が求めることができる。 この方程式の根は、二つあるので、t の値の侯補も、二つあることになるが、 その内 \thita の絶対値の小さい方を取るので、 結局、 t = \frac{\sign{\phi}}{|\phi|+\sqrt{\phi^2+1}} 但し、 \phi = \frac{a_qq-a_pp}{2a_pq} が選ばれる。 これより、 \cos{\thita} = \frac{1}{\sqrt{1+t^2}} \sin{\thita} = t\cos{\thita} を得ることができる。 # 実際には、\thita を計算することはない ( \thita は数学的に利用するだけ ) # 実際に利用するのは、\sin{\thita} と \cos{\thita} だけ # \thita が嫌な理由 # => \arctan が「(数値的に..)嫌」な関数だから !! == [課題] ヤコビ法は、複雑なので、二回 ( 今日 : 09/27 と、来々週 : 10/11 の二日 ) に分けて説明する。 # したがって、課題 (提出) も、一つ 課題番号は、10/11 締切は、10/17 プログラムは、text p.275 例題は、次回、説明する