|
DGESVD/ZGESVD (CLAPACK) |
DGESVD/ZGESVD
一般長方行列の特異値分解をして、特異値ベクトルを計算する。 A = U * Σ * VH
＜関数＞
int dgesvd_(char *jobu, char *jobvt, integer *m, integer *n, doublereal *a,
integer *lda, doublereal *s, doublereal *u, integer * ldu, doublereal *vt,
integer *ldvt, doublereal *work, integer *lwork, integer *info);
＜引数＞
char jobu - (input) 行列[U]の列(左特異値ベクトル）を計算し格納するかどうかを指示する。
指定できる選択肢は以下の通りである。（以下に載っていない値は規定外である）
='A': 行列[U]のm個のすべての列(左特異値ベクトル)を配列Uに格納する。
='N': 行列[U]の列(左特異値ベクトル)は計算しない。
='O':
行列[U]の最初のmin(m,n)個の列(左特異値ベクトル)で配列Aを上書きする。
すなわち、特異値に対応する左特異値ベクトルを配列Aに上書き。
='S':
行列[U]の最初のmin(m,n)個の列(左特異値ベクトル)だけが配列Uに格納される。
すなわち、特異値に対応する左特異値ベクトルを配列Uに格納する。
char jobvt - (input) 行列[V]Hの行(右特異値ベクトル）を計算し格納するかどうかを指示する。
='A': 行列[V]Hのn個のすべての行(右特異値ベクトル)を配列VTに格納する。
='N': 行列[V]Hの列(右特異値ベクトル)は計算しない。
='O':
行列[V]Hの最初のmin(m,n)個の行(右特異値ベクトル)で配列Aを上書きする。
すなわち、特異値に対応する右特異値ベクトルを配列Aに上書き。
='S':
行列[V]Hの最初のmin(m,n)個の行(右特異値ベクトル)だけが配列VTに格納される。
すなわち、特異値に対応する右特異値ベクトルを配列VTに格納する。
long int m - (input)入力する行列[A]の行数。M≧0。
long int n - (input)入力する行列[A]の列数。N≧0。
double a - (input/output)配列形式はa[lda,n]。
入力時: (m,n)の行列。
出力時:
jobu = 'O'もしくは'o'ならば、行列[U]の最初のmin(m,n)個の列(左特異値ベクトル)で配列Aは上書きされる。
jobvr = 'O'もしくは'o'ならば、行列[V]Hの最初のmin(m,n)個の行(右特異値ベクトル)で配列Aは上書きされる。
jobuもjobvrも'O'もしくは'o'ならば、Aの内容は上書きされる。
long int lda - (input)配列Aの第一次元(のメモリ格納数)。lda≧max(1,m) つまり、1〜mの中で一番大きい数よりも、ldaが大きくなければならない。
通常はlda=mで良い。
double s - (output) min(m,n)個の行列Aの特異値。配列の次元数はmin(m,n)。降順で並べられている。特異値は対角行列Σの対角要素である。
double u - (output) 配列uの内容はjobuに指定した値によって変わる。詳細は以下の通りである。
jobu='A' or 'a': 行列[U]のサイズはm×m。直行行列かユニタリ行列Uが配列Uに入る。
jobu='N' or 'n': 行列[U]のサイズは0。行列[U]は使われない。
jobu='O or 'o'': 行列[U]のサイズは0。行列[U]は使われない。
jobu='S' or 's': 行列[U]のサイズはm×min(m,n)。直行行列かユニタリ行列Uの最初のmin(m,n)個の列(左特異ベクトル)が行列[U]に入る。
long int ldu - (input) 配列uの第一次元(のメモリ格納数)。ldu≧1。jobu = 'A', 'a' or 'S' or 's'のとき、ldu≧max(1,m)。
double vt - (output) 配列の大きさは、ldvt×nの大きさで取ること。配列vtの内容はjobvtに指定した値によって変わる。詳細は以下の通りである。
jobvt='A' or 'a': 行列[VT]のサイズはn×n。配列vtには直行行列かユニタリ行列VTが入る。
jobvt='N' or 'n': 行列[VT]のサイズは0。
jobvt='O' or 'o': 行列[VT]のサイズは0。
jobvt='S' or 's': 直行行列かユニタリ行列VTの最初のmin(m,n)個の行(右特異ベクトル)が行列[VT]に入る。
long int ldvt - (input) 配列vtの第一次元(のメモリ格納数)。詳細は以下に従う。
jobvt = 'A' or 'a' のとき、ldvt≧max(1,n)
jobvt = 'S' or 's' のとき、ldvt≧min(m,n) かつ ldvt≧1。
jobvtが他の値のとき、ldvt≧1。
double work - (working/output) 次元数lworkの計算用配列。
long int lwork - (input) 配列workの次元数。lwork≧1であり、以下の式に従う。
実数関数(DGESVD)の場合:
lwork≧max(3*min(m,n)+max(m,n), 5*min(m,n)-4)。
しかし、lwork≧max(3*min(m,n)+max(m,n),
5*min(m,n))にしないと動かない事を確認。多めにとった方が良さそうです。
複素数関数(ZGESVD)の場合: lwork≧2*min(m,n)+max(m,n)
long int info - (output)
info = 0: 正常に終了した事を示す。
info < 0: info = -i ならば、i番目の引数の値が間違えていることを示す。
info > 0: 計算が収束しなかかった事を示す。
を特異値分解する。
|
#include <stdio.h> #define ROW 4 #define COL 4 double A[ROW*COL]; int main(void) { static long int i; char jobu = 'N'; char jobvt = 'N'; static long int m = ROW, n = COL; static long int lda=ROW, ldu=1, ldvt=1, lwork=20, info; static double s[ROW], u[1], vt[1], work[20]; A[0]=0.5;A[1]=0.5;A[2]=0.5;A[3]=0.5; A[4]=0.5;A[5]=0.5;A[6]=-0.5;A[7]=-0.5; A[8]=0.5;A[9]=-0.5;A[10]=0.5;A[11]=-0.5; A[12]=0.5;A[13]=-0.5;A[14]=-0.5;A[15]=0.5; dgesvd_(&jobu, &jobvt, &m, &n, A, &lda, s, u, &ldu, vt, &ldvt, work, &lwork, &info); for(i=0;i<ROW;++i) printf("%lf\n",s[i]); return 0; } |
Copyright (C) 2001 Keisuke ABE. All Rights Reserved.