Zgejsv Zgesdd Zgesvd Zgesvdq Zgesvdx Zggsvd3

Zgesvdx() Sub Zgesvdx ( Jobu As  String, Jobvt As  String, Range As  String, M As  Long, N As  Long, A() As  Complex, Vl As  Double, Vu As  Double, Il As  Long, Iu As  Long, Ns As  Long, S() As  Double, U() As  Complex, Vt() As  Complex, Info As  Long  ) (エキスパートドライバ) 特異値分解 (SVD) (複素行列) 目的 本ルーチンはm×n複素行列Aの特異値分解(SVD), および, 必要により左および／または右特異ベクトルを求める. SVDは次のように表される. A = U * SIGMA * V^H ここで, SIGMAはmin(m, n)個の対角要素を除き0のm×n行列, Uはm×mユニタリ行列, Vはn×nユニタリ行列である. SIGMAの対角要素がAの特異値である. 特異値は非負の実数で, 降順に返される. UおよびVのはじめのmin(m, n)列はAの左および右特異ベクトルである. Zgesvdxは, SVDを求めるのに固有値問題を使用する. これにより, 特異値およびベクトルの一部だけを計算することができる. 本ルーチンではVではなくV^Hを返すので注意せよ. 引数 [in] Jobu 行列Uの全部あるいは一部を計算するか指定する. = "V": U(左特異ベクトル)のはじめのmin(M, N)列, または, Rangeで指定された一部の特異値に対応する列を配列U()に返す. = "N": Uの列(左特異ベクトル)を計算しない. [in] Jobvt 行列V^Hの全部あるいは一部を計算するか指定する. = "V": V^H(右特異ベクトル)のはじめのmin(M, N)行, または, Rangeで指定された一部の特異値に対応する行を配列Vt()に返す. = "N": V^Hの行(右特異ベクトル)を計算しない. [in] Range = "A": すべての特異値を求める. = "V": 半開区間 (Vl, Vu] のすべての特異値を求める. = "I": Il番目からIu番目までの特異値を求める. [in] M 入力行列Aの行数. (M >= 0) (M = 0 の場合, 処理を行わずに戻る) [in] N 入力行列Aの列数. (N >= 0) (N = 0 の場合, 処理を行わずに戻る) [in,out] A() 配列 A(LA1 - 1, LA2 - 1) (LA1 >= M, LA2 >= N) [in] M×N 行列 A. [out] A()の内容は壊される. [in] Vl Range = "V": 求める特異値の範囲の下限. (Vu > Vl) Range = "A" または "I": 参照されない. [in] Vu Range = "V": 求める特異値の範囲の上限. (Vu > Vl) Range = "A" または "I": 参照されない. [in] Il Range = "I"] 求める最小特異値のインデックス. (1 <= Il <= Iu <= min(M,N), ただし min(M,N) > 0) Range = "A" または "V": 参照されない. [in] Iu Range = "I": 求める最大特異値のインデックス. (1 <= Il <= Iu <= min(M,N), ただし min(M,N) > 0) Range = "A" または "V": 参照されない. [out] Ns 求められた特異値の総数. (0 <= Ns <= min(M,N)) Range = "A": Ns = min(M, N). Range = "I": Ns = Iu - Il + 1. [out] S() 配列 S(LS - 1) (LS >= min(M, N)) Aの特異値 (S(i) >= S(i+1) となるように並べ替えられる). [out] U() 配列 U(LU1 - 1, LU2 - 1) (LU1 >= M, LU2 >= NS) Jobu = "V": Rangeで指定された特異値に対応するU(左特異ベクトル, 列ごとに格納)の列がU()に入る. Jobu = "N": U()は参照されない. 注: LU2 >= NS となるようにすること. Range = "V"の場合, NSの正確な値は前もってわからないが上限値を使用しなければならない. [out] Vt() 配列 Vt(LVt1 - 1, LVt2 - 1) (LVt1 >= NS, LVt2 >= N) Jobvt = "V": Rangeで指定された特異値に対応するV^H(右特異ベクトル, 行方向に格納)の行がVt()に入る. Jobvt = "N": Vt()は参照されない. 注: LVt1 >= NS となるようにすること. Range = "V"の場合, NSの正確な値は前もってわからないが上限値を使用しなければならない. [out] Info = 0: 正常終了. = -1: パラメータ Jobu の誤り. (Jobu <> "A", "S", "O"および"N") = -2: パラメータ Jobvt の誤り. (Jobvt <> "A", "S", "O"および"N") = -3: パラメータ Range の誤り. (Range != "A", "V"および"I") = -4: パラメータ M の誤り. (M < 0) = -5: パラメータ N の誤り. (N < 0) = -6: パラメータ A() の誤り. = -7: パラメータ Vl の誤り. (Vl < 0) = -8: パラメータ Vu の誤り. (Vu < Vl) = -9: パラメータ Il の誤り. (Il < 1 または Il > max(1, min(M, N))) = -10: パラメータ Iu の誤り. (Iu < min(min(M, N), Il) または Iu > min(M, N)) = -12: パラメータ S() の誤り. = -13: パラメータ U() の誤り. = -14: パラメータ Vt() の誤り. = i > 0: Dbdsvdx/Dstevxにおいてi個の固有ベクトルが収束しなかった. = 2*N + 1: Dbdsvdxにおいて内部エラーが起きた. 出典 LAPACK 使用例 行列Aの特異値, 左および右特異ベクトルを求める. ただし, ( 0.20-0.11i -0.93-0.32i 0.81+0.37i ) A = ( -0.80-0.92i -0.29+0.86i 0.64+0.51i ) ( 0.71+0.59i -0.15+0.19i 0.20+0.94i ) とする. Sub Ex_Zgesvdx() Const M = 3, N = 3 Dim A(M - 1, N - 1) As Complex, U(M - 1, N - 1) As Complex, Vt(N - 1, N - 1) As Complex Dim S(N - 1) As Double, Info As Long Dim Vl As Double, Vu As Double, Il As Long, Iu As Long, Ns As Long A(0, 0) = Cmplx(0.2, -0.11): A(0, 1) = Cmplx(-0.93, -0.32): A(0, 2) = Cmplx(0.81, 0.37) A(1, 0) = Cmplx(-0.8, -0.92): A(1, 1) = Cmplx(-0.29, 0.86): A(1, 2) = Cmplx(0.64, 0.51) A(2, 0) = Cmplx(0.71, 0.59): A(2, 1) = Cmplx(-0.15, 0.19): A(2, 2) = Cmplx(0.2, 0.94) Call Zgesvdx("V", "V", "A", M, N, A(), Vl, Vu, Il, Iu, Ns, S(), U(), Vt(), Info) Debug.Print "Singular values =", S(0), S(1), S(2) Debug.Print "U =" Debug.Print Creal(U(0, 0)), Cimag(U(0, 0)), Creal(U(0, 1)), Cimag(U(0, 1)) Debug.Print Creal(U(1, 0)), Cimag(U(1, 0)), Creal(U(1, 1)), Cimag(U(1, 1)) Debug.Print Creal(U(2, 0)), Cimag(U(2, 0)), Creal(U(2, 1)), Cimag(U(2, 1)) Debug.Print Creal(U(0, 2)), Cimag(U(0, 2)) Debug.Print Creal(U(1, 2)), Cimag(U(1, 2)) Debug.Print Creal(U(2, 2)), Cimag(U(2, 2)) Debug.Print "V^H =" Debug.Print Creal(Vt(0, 0)), Cimag(Vt(0, 0)), Creal(Vt(0, 1)), Cimag(Vt(0, 1)) Debug.Print Creal(Vt(1, 0)), Cimag(Vt(1, 0)), Creal(Vt(1, 1)), Cimag(Vt(1, 1)) Debug.Print Creal(Vt(2, 0)), Cimag(Vt(2, 0)), Creal(Vt(2, 1)), Cimag(Vt(2, 1)) Debug.Print Creal(Vt(0, 2)), Cimag(Vt(0, 2)) Debug.Print Creal(Vt(1, 2)), Cimag(Vt(1, 2)) Debug.Print Creal(Vt(2, 2)), Cimag(Vt(2, 2)) Debug.Print "Ns =", Ns, "Info =", Info End Sub 実行結果 Singular values = 2.07084030821889 1.23513084760163 0.901483337149816 U = -0.322093837376809 0.402732043434636 9.95958553097354E-02 0.309610800223195 0.25572227156877 0.689932737910445 -0.63162108895219 3.70329357348017E-02 -0.438560869464762 1.80488745428176E-02 -7.00863751874233E-02 0.699280401316292 0.711819740685389 -0.348707521966052 -9.32902099996817E-02 0.222663951646943 -0.551156386049191 -9.82856061041459E-02 V^H = -0.603023132575377 0 0.36655126147564 0.394522570687146 0.663818465313821 0 0.134961467338043 -0.149629139353428 -0.44238913491109 0 -0.297134277147366 -0.762299060909463 0.16075433957414 -0.566138903287951 0.366911928609536 -0.61977189734833 0.331437452593181 -0.158279598849071 Ns = 3 Info = 0