Dgges Dgges_r Dggesx Dggesx_r Dggev Dggevx

Dggevx() Sub Dggevx ( Balanc As  String, Jobvl As  String, Jobvr As  String, Sense As  String, N As  Long, A() As  Double, B() As  Double, Alphar() As  Double, Alphai() As  Double, Beta() As  Double, Vl() As  Double, Vr() As  Double, Ilo As  Long, Ihi As  Long, LScale() As  Double, RScale() As  Double, AbNrm As  Double, BbNrm As  Double, RConde() As  Double, RCondv() As  Double, Info As  Long  ) (エキスパートドライバ) 一般化固有値問題 (一般行列) 目的 本ルーチンはn×n一般実行列のペア(A, B)の一般化固有値, および, 必要により左および／または右一般化固有ベクトルを求める. また, 必要により固有値および固有ベクトルの条件を改善するための均衡化の変換(Ilo, Ihi, Lscale, Rscale, AbNrmおよびBbNrm), 固有値の条件数の逆数(RConde), および, 右固有ベクトルの条件数の逆数(RCondv)を求める. 行列のペア(A, B)の一般化固有値は, A - λB が特異となるようなスカラーλあるいは比 α/β=λである. β = 0 あるいは両方共0の場合に適当な解釈ができるため, これは通常ペア(α, β)で表される. (A, B)の一般化固有値λ(j)に対応する右一般化固有ベクトルv(j)は次式を満たす. A * v(j) = λ(j) * B * v(j) (A, B)の一般化固有値λ(j)に対応する左一般化固有ベクトルu(j)は次式を満たす. u(j)^H * A = λ(j) * u(j)^H * B ここで, u(j)^Hはu(j)の共役転置である. 引数 [in] Balanc 行うべき均衡化の方法を指定する. = "N": 対角スケーリングおよび入替を行わない. = "P": 入替えのみを行う. = "S": 対角スケーリングのみを行う. = "B": 対角スケーリングおよび入替を行う. 条件数の逆数は入替および／またはスケーリング後の行列について求められる. 入替により条件数は変わらない(丸め誤差がなければ)が, スケーリングにより変わる. [in] Jobvl = "N": 左一般化固有ベクトルを求めない. = "V": 左一般化固有ベクトルを求める. [in] Jobvr = "N": 右一般化固有ベクトルを求めない. = "V": 右一般化固有ベクトルを求める. [in] Sense どの条件数の逆数を計算するかを決定する. = "N": 計算しない. = "E": 固有値についてのみ計算する. = "V": 固有ベクトルについてのみ計算する. = "B": 固有値および固有ベクトルについて計算する. [in] N 行列A, B, VLおよびVRの行および列数. (N >= 0) (N = 0 の場合, 処理を行わずに戻る) [in,out] A() 配列 A(LA1 - 1, LA2 - 1) (LA1 >= N, LA2 >= N) [in] ペア(A, B)の行列 A. [out] A()は上書きされる. Jobvl = "V"あるいはJobvr = "V"あるいは両方の場合, 均衡化された入力AおよびBの実シュール形の最初の部分が入る. [in,out] B() 配列 B(LB1 - 1, LB2 - 1) (LB1 >= N, LB2 >= N) [in] ペア(A, B)の行列 B. [out] B()は上書きされる. Jobvl = "V"あるいはJobvr = "V"あるいは両方の場合, 均衡化された入力AおよびBの実シュール形の続きの部分が入る. [out] Alphar() 配列 Alphar(LAlphar - 1) (LAlphar >= N) [out] Alphai() 配列 Alphai(LAlphai - 1) (LAlphai >= N) [out] Beta() 配列 Beta(LBeta - 1) (LBeta >= N) (Alphar(j) + Alphai(j)*i)/Beta(j) (j = 0〜N-1) は一般化固有値である. Alphai(j)が0の場合, j番目の固有値は実数である. 正数の場合, j番目とj+1番目の固有値は複素共役対であり, Alphai(j+1)は負である. 注: 商Alphar(j)/Beta(j)およびAlphai(j)/Beta(j)はオーバーフローあるいはアンダーフローし易く, Beta(j)が0になることさえある. 従って, 比α/βを単純に計算することは避けなければならない. しかし, AlpharとAlphaiは常にnorm(A)より小さく, 通常は同程度の大きさである. また, Betaは常にnorm(B)より小さく, 通常は同程度の大きさである. [out] Vl() 配列 Vl(LVl1 - 1, LVl2 - 1) (LVl1 >= N, LVl2 >= N) Jobvl = "V": 左固有ベクトルu(j)が固有値と同じ順にVl()の列に入る.   j番目の固有値が実数の場合, u(j) = (Vl()のj番目の列) である.   j番目とj+1番目の固有値が複素共役対の場合, u(j) = (Vl()のj番目の列) + (Vl()のj+1番目の列)*i, u(j+1) = (Vl()のj番目の列) - (Vl()のj+1番目の列)*i である.   固有ベクトルは, 最大要素が |実数部| + |虚数部| = 1 となるように正規化される. Jobvl = "N": Vl()は参照されない. [out] Vr() 配列 Vr(LVr1 - 1, LVr2 - 1) (LVr1 >= N, LVr2 >= N) Jobvr = "V": 右固有ベクトルv(j)が固有値と同じ順にVr()の列に入る.   j番目の固有値が実数の場合, v(j) = (Vr()のj番目の列) である.   j番目とj+1番目の固有値が複素共役対の場合, v(j) = (Vr()のj番目の列) + (Vr()のj+1番目の列)*i, v(j+1) = (Vr()のj番目の列) - (Vr()のj+1番目の列)*i である.   固有ベクトルは, 最大要素が |実数部| + |虚数部| = 1 となるように正規化される. Jobvr = "N": Vr()は参照されない. [out] Ilo [out] Ihi IloとIhiは, i > j かつ j = 0〜Ilo-2 または i = Ihi〜N-1 に対して, A(i, j) = 0, B(i, j) = 0 となるような整数値である. balanc = "N"あるいは"S"の場合, Ilo = 1, Ihi = N となる. [out] LScale() 配列 LScale(LLScale - 1) (LLScale >= N) 左からAとBに適用される入替とスケーリング因子の詳細. Pl(j)を行jと入替を行った行のインデックス, Dl(j)を行jに適用したスケーリング因子としたとき, 次のようになる.   LScale(j) = Pl(j) (j = 0 〜 Ilo-2)     = Dl(j) (j = Ilo-1 〜 Ihi-1)     = Pl(j) (j = Ihi 〜 N-1) 入替の順番はN-1からIhi, そして0からIlo-2である. [out] RScale() 配列 RScale(LRScale - 1) (LRScale >= N) 右からAとBに適用される入替とスケーリング因子の詳細. Pr(j)を列jと入替を行った列のインデックス, Dr(j)を列jに適用したスケーリング因子としたとき, 次のようになる.   RScale(j) = Pr(j) (j = 0 〜 Ilo-2)     = Dr(j) (j = Ilo-1 〜 Ihi-1)     = Pr(j) (j = Ihi 〜 N-1) 入替の順番はN-1からIhi, そして0からIlo-2である. [out] AbNrm 均衡化された行列Aの1ノルム. [out] BbNrm 均衡化された行列Bの1ノルム. [out] RConde() 配列 RConde(LRConde - 1) (LRConde >= N) Sense = "E"または"B": 固有値の条件数の逆数が配列の要素に順に格納される. 複素共役対の固有値については, RConde()の連続する2つの要素に同一の値が設定される. 従って, RConde(j), RCondv(j), および, Vl()とVr()のj列はすべてj番目の固有値対に対応する. Sense = "N"または"V": RConde()は参照されない. [out] RCondv() 配列 RCondv(LRCondv - 1) (LRCondv >= N) Sense = "V"または"B": 固有ベクトルの条件数の逆数の推定値が配列の要素に順に格納される. 複素固有ベクトルについては, RCondv()の連続する2つの要素に同一の値が設定される. RCondv(j)を求めるために固有値を並べ替えることができない場合, RCondv(j)は0に設定される. これは, 真の値が非常に小さい場合にのみ起こりうる. Sense = "N"または"E": RCondv()は参照されない. [out] Info = 0: 正常終了. = -1: パラメータ Balanc の誤り. (Balanc <> "N", "P", "S"および"B") = -2: パラメータ Jobvl の誤り. (Jobvl <> "V"および"N") = -3: パラメータ Jobvr の誤り. (Jobvr <> "V"および"N") = -4: パラメータ Sense の誤り. (Sense <> "N", "E", "V"および"B") = -5: パラメータ N の誤り. (N < 0) = -6: パラメータ A() の誤り. = -7: パラメータ B() の誤り. = -8: パラメータ Alphar() の誤り. = -9: パラメータ Alphai() の誤り. = -10: パラメータ Beta() の誤り. = -11: パラメータ Vl() の誤り. = -12: パラメータ Vr() の誤り. = -15: パラメータ LScale() の誤り. = -16: パラメータ RScale() の誤り. = -19: パラメータ RConde() の誤り. = -20: パラメータ RCondv() の誤り. = i (0 < i <= N): QZ反復が収束しなかった. 固有ベクトルは計算されないが, Alphar(j), Alphai(j), Beta(j) (j = i〜N-1) は正しい値である. = N+1: DhgeqzにおいてQZ反復の収束以外のエラーが起きた. = N+2: Dtgevcにおいてエラーがあった. 出典 LAPACK 使用例 行列のペア(A, B)の一般化固有値および左および右一般化固有ベクトルを求める. ただし, ( 0.20 -0.11 -0.93 ) ( -0.58 -0.79 0.82 ) A = ( -0.32 0.81 0.37 ), B = ( 0.77 0.71 -0.55 ) ( -0.80 -0.92 -0.29 ) ( -1.36 -1.22 1.66 ) とする. Sub Ex_Dggevx() Const N = 3 Dim A(N - 1, N - 1) As Double, B(N - 1, N - 1) As Double Dim Alphar(N - 1) As Double, Alphai(N - 1) As Double, Beta(N - 1) As Double Dim Vl(N - 1, N - 1) As Double, Vr(N - 1, N - 1) As Double, Info As Long Dim Ilo As Long, Ihi As Long, LScale(N - 1) As Double, RScale(N - 1) As Double Dim AbNrm As Double, BbNrm As Double Dim RConde(N - 1) As Double, RCondv(N - 1) As Double A(0, 0) = 0.2: A(0, 1) = -0.11: A(0, 2) = -0.93 A(1, 0) = -0.32: A(1, 1) = 0.81: A(1, 2) = 0.37 A(2, 0) = -0.8: A(2, 1) = -0.92: A(2, 2) = -0.29 B(0, 0) = -0.58: B(0, 1) = -0.79: B(0, 2) = 0.82 B(1, 0) = 0.77: B(1, 1) = 0.71: B(1, 2) = -0.55 B(2, 0) = -1.36: B(2, 1) = -1.22: B(2, 2) = 1.66 Call Dggevx("N", "V", "V", "N", N, A(), B(), Alphar(), Alphai(), Beta(), Vl(), Vr(), Ilo, Ihi, LScale(), RScale(), AbNrm, BbNrm, RConde(), RCondv(), Info) Debug.Print "Eigenvalues =" Debug.Print " (r)", Alphar(0) / Beta(0), Alphar(1) / Beta(1), Alphar(2) / Beta(2) Debug.Print " (i)", Alphai(0) / Beta(0), Alphai(1) / Beta(1), Alphai(2) / Beta(2) Debug.Print "Eigenvectors (L) =" Debug.Print Vl(0, 0), Vl(0, 1), Vl(0, 2) Debug.Print Vl(1, 0), Vl(1, 1), Vl(1, 2) Debug.Print Vl(2, 0), Vl(2, 1), Vl(2, 2) Debug.Print "Eigenvectors (R) =" Debug.Print Vr(0, 0), Vr(0, 1), Vr(0, 2) Debug.Print Vr(1, 0), Vr(1, 1), Vr(1, 2) Debug.Print Vr(2, 0), Vr(2, 1), Vr(2, 2) Debug.Print "Info =", Info End Sub 実行結果 Eigenvalues = (r) 1.10765548065266 1.10765548065266 -1.74329621491312 (i) 1.40252005479606 -1.40252005479606 0 Eigenvectors (L) = 0.828456655019124 0.171543344980876 -1 -0.174134648086565 0.261212076890975 -0.66437599243234 -0.601718475294568 -0.175396675418769 4.15956472218708E-02 Eigenvectors (R) = -9.54330235945138E-02 0.502674461662748 -1 -0.945775638065751 5.42243619342491E-02 0.224694817103864 -0.659251386021417 5.30923334667921E-02 -0.954838851339161 Info = 0