Dgejsv Dgesdd Dgesvd Dgesvdq Dgesvdx Dggsvd3

Dgesvd() Sub Dgesvd ( Jobu As  String, Jobvt As  String, M As  Long, N As  Long, A() As  Double, S() As  Double, U() As  Double, Vt() As  Double, Info As  Long  ) (シンプルドライバ) 特異値分解 (SVD) 目的 本ルーチンはm×n実行列Aの特異値分解(SVD), および, 必要により左および／または右特異ベクトルを求める. SVDは次のように表される. A = U * Σ * V^T ここで, Σはmin(m, n)個の対角要素を除き0のm×n行列, Uはm×m直交行列, Vはn×n直交行列である. Σの対角要素がAの特異値である. 特異値は非負の実数で, 降順に返される. UおよびVのはじめのmin(m, n)列はAの左および右特異ベクトルである. 本ルーチンではVではなくV^Tを返すので注意せよ. 引数 [in] Jobu 行列Uの全部あるいは一部を計算するか指定する. = "A": UのM列全てを配列U()に返す. = "S": Uのはじめのmin(M, N)列(左特異ベクトル)を配列U()に返す. = "O": Uのはじめのmin(M, N)列(左特異ベクトル)を配列A()に上書きする. = "N": Uの列(左特異ベクトル)を計算しない. 注 - JobvtとJobuを両方"O"にはできない. [in] Jobvt 行列V^Tの全部あるいは一部を計算するか指定する. = "A": V^TのN行全てを配列Vt()に返す. = "S": V^Tのはじめのmin(M, N)行(右特異ベクトル)を配列Vt()に返す. = "O": V^Tのはじめのmin(M, N)行(右特異ベクトル)を配列A()に上書きする. = "N": V^Tの行(右特異ベクトル)を計算しない. 注 - JobvtとJobuを両方"O"にはできない. [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] Jobu = "O": A()はU(左特異ベクトル, 列ごとに格納)のはじめのmin(M, N)列で上書きされる.   Jobvt = "O": A()はV^T(右特異ベクトル, 行ごとに格納)のはじめのmin(M, N)行で上書きされる.   Jobu <> "O" かつ Jobvt <> "O": A()の内容は壊される. [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 >= M (Jobu = "A") または LU2 >= min(M, N) (Jobu = "S")) Jobu = "A": U()にM×M直交行列Uが入る. Jobu = "S": U()にU(左特異ベクトル, 列ごとに格納)のはじめのmin(M, N)列が入る. Jobu = "N"あるいは"O": U()は参照されない. [out] Vt() 配列 Vt(LVt1 - 1, LVt2 - 1) (LVt1 >= N (JobVt = "A") または LVt1 >= min(M, N) (Jobvt = "S"), LVt2 >= N) Jobvt = "A": Vt()にN×N直交行列V^Tが入る. Jobvt = "S": Vt()にV^T(右特異ベクトル, 行ごとに格納)のはじめのmin(M, N)行が入る. Jobvt = "N"あるいは"O": Vt()は参照されない. [out] Info = 0: 正常終了. = -1: パラメータ Jobu の誤り. (Jobu <> "A", "S", "O"および"N") = -2: パラメータ Jobvt の誤り. (Jobvt <> "A", "S", "O"および"N") = -3: パラメータ M の誤り. (M < 0) = -4: パラメータ N の誤り. (N < 0) = -5: パラメータ A() の誤り. = -6: パラメータ S() の誤り. = -7: パラメータ U() の誤り. = -8: パラメータ Vt() の誤り. > 0: Dbdsqrで収束しなかった場合, Infoは中間結果の二重対角形Bの上副対角要素のうち何個が0に収束しなかったを示す. 出典 LAPACK 使用例 行列Aの特異値, 左および右特異ベクトルを求める. ただし, ( 1 6 11 ) ( 2 7 12 ) A = ( 3 8 13 ) ( 4 9 14 ) ( 5 10 15 ) である. Sub Ex_Dgesvd() Const M = 5, N = 3 Dim A(M - 1, N - 1) As Double, U(M - 1, N - 1) As Double, Vt(N - 1, N - 1) As Double Dim S(N - 1) As Double, Info As Long A(0, 0) = 1: A(0, 1) = 6: A(0, 2) = 11 A(1, 0) = 2: A(1, 1) = 7: A(1, 2) = 12 A(2, 0) = 3: A(2, 1) = 8: A(2, 2) = 13 A(3, 0) = 4: A(3, 1) = 9: A(3, 2) = 14 A(4, 0) = 5: A(4, 1) = 10: A(4, 2) = 15 Call Dgesvd("S", "A", M, N, A(), S(), U(), Vt(), Info) Debug.Print "Singular values =" Debug.Print S(0), S(1), S(2) Debug.Print "U =" Debug.Print U(0, 0), U(0, 1), U(0, 2) Debug.Print U(1, 0), U(1, 1), U(1, 2) Debug.Print U(2, 0), U(2, 1), U(2, 2) Debug.Print U(3, 0), U(3, 1), U(3, 2) Debug.Print U(4, 0), U(4, 1), U(4, 2) Debug.Print "V^T =" Debug.Print Vt(0, 0), Vt(0, 1), Vt(0, 2) Debug.Print Vt(1, 0), Vt(1, 1), Vt(1, 2) Debug.Print Vt(2, 0), Vt(2, 1), Vt(2, 2) Debug.Print "Info =", Info End Sub 実行結果 Singular values = 35.1272233335747 2.46539669691652 1.26909062514315E-15 U = -0.354557057037681 0.688686643768252 0.570044169987079 -0.398696369998832 0.375554529395871 -0.74547800114907 -0.442835682959984 6.24224150234906E-02 -0.170157929158505 -0.486974995921135 -0.25070969934889 0.296573181815902 -0.531114308882286 -0.56384181372127 4.90185785045937E-02 V^T = -0.201664911192694 -0.516830501392304 -0.831996091591915 -0.890317132783019 -0.257331626824052 0.375653879134918 0.408248290463864 -0.816496580927726 0.408248290463863 Info = 0