?gesdd

Compute the SVD of a rectangular matrix, and optionally the left and right singular vectors, using a divide-and-conquer method.

Interface Definition

C interface:

void sgesdd_(const char *jobz, const int *m, const int *n, float *a, const int *lda, float *s, float *u, const int *ldu, float *vt, const int *ldvt, float *work, const int *lwork, int *iwork, int *info);

void dgesdd_(const char *jobz, const int *m, const int *n, double *a, const int *lda, double *s, double *u, const int *ldu, double *vt, const int *ldvt, double *work, const int *lwork, int *iwork, int *info);

void cgesdd_(const char *jobz, const int *m, const int *n, float _Complex *a, const int *lda, float *s, float _Complex *u, const int *ldu, float _Complex *vt, const int *ldvt, float _Complex *work, const int *lwork, float *rwork, int *iwork, int *info);

void zgesdd_(const char *jobz, const int *m, const int *n, double _Complex *a, const int *lda, double *s, double _Complex *u, const int *ldu, double _Complex *vt, const int *ldvt, double _Complex *work, const int *lwork, double *rwork, int *iwork, int *info);

Fortran interface:

SGESDD(JOBZ, M, N, A, LDA, S, U, LDU, VT, LDVT, WORK, LWORK, IWORK, INFO)

DGESDD(JOBZ, M, N, A, LDA, S, U, LDU, VT, LDVT, WORK, LWORK, IWORK, INFO)

CGESDD(JOBZ, M, N, A, LDA, S, U, LDU, VT, LDVT, WORK, LWORK, RWORK, IWORK, INFO)

ZGESDD(JOBZ, M, N, A, LDA, S, U, LDU, VT, LDVT, WORK, LWORK, RWORK, IWORK, INFO)

Parameters

Parameter	Type	Description	Input/Output
jobz	Character	'N': computes only eigenvalues and does not return U or Vt. 'S': returns the first min(M, N) columns of U and the first min(M, N) rows of Vt. 'O': If M ≥ N, the first N columns of U are overwritten on matrix A, and all rows of Vt are returned in array Vt. 'A': All M columns of U and all N rows of Vt are returned in arrays U and Vt, respectively.	Input
m	Integer	Dimension of the matrix. m ≥ 0.	Input
n	Integer	Matrix dimension. n ≥ 0.	Input
a	A single-precision floating-point array for sgesdd A double-precision floating-point array for dgesdd A single-precision complex array for cgesdd A double-precision complex array for zgesdd	An array of size *MN. On entry: a matrix of size MN. On exit: If jobz* = 'O' and M ≥ N, A is overwritten by the first N columns of U (left singular vectors, column-major order). If jobz != 'O', A is destroyed.	Input/Output
lda	Integer	Leading dimension of matrix A. lda ≥ max(1, M).	Input
s	A single-precision floating-point array for sstedc A double-precision floating-point array for dstedc A single-precision floating-point array for cstedc A double-precision floating-point array for zstedc	The size is min(M, N). Singular values of matrix A, sorted in descending order.	Output
u	A single-precision floating-point array for sgesdd A double-precision floating-point array for dgesdd A single-precision complex array for cgesdd A double-precision complex array for zgesdd	A two-dimensional array with a size of *LDUUCOL. If jobz = 'A' or 'O' and M < N, UCOL = M. u contains the M x M orthogonal matrix. If jobz = 'S', UCOL = min(M, N). uz contains the first min(M, N) columns of the left singular vectors U in the column-major order. If jobz = 'O' and M ≥ N, or jobz = 'N'**, this parameter is ignored.	Output
ldu	Integer	Leading dimension of matrix u. ldu ≥ 1. If jobz = 'S' or 'A', or jobz = 'O' and M < N, LDU ≥ M.	Input
vt	A single-precision floating-point array for sgesdd A double-precision floating-point array for dgesdd A single-precision complex array for cgesdd A double-precision complex array for zgesdd	A two-dimensional array with a size of *LDVTN. If jobz = 'A' or 'O' and M < N, UCOL = M and it contains the M x M orthogonal matrix. If jobz = 'S', it contains the first min(M, N) rows of the right singular vectors in the row-major order. If jobz = 'O' and M ≥ N, or jobz = 'N'**, this parameter is ignored.	Output
ldvt	Integer	Leading dimension of the matrix vt. ldvt ≥ 1. If jobz = 'A' or 'O' and M < N, ldvt ≥ N. If jobz = 'S', ldvt ≥ min (M, N).	Input
work	A single-precision floating-point array for sgesdd A double-precision floating-point array for dgesdd A single-precision complex array for cgesdd A double-precision complex array for zgesdd	Array of dimension max(1, LWORK). On exit: If info = 0, the first element returns the optimal lwork.	Output
lwork	Integer	Dimension of matrix work. lwork ≥ 1. If lwork = 1, the query is performed by default. The optimal size is returned in the first element of work, and no other work parameters are checked. If JOBZ = 'N', LWORK ≥ 3mn + max(mx, 7mn). If JOBZ = 'O', LWORK ≥ 3mn + max(mx, 5mnmn + 4mn). If JOBZ = 'S', LWORK ≥ 4mnmn + 7mn. If JOBZ = 'A', lwork* ≥ 4mnmn + 6*mn + mx.	Input
rwork (only for complex types)	A single-precision floating-point array for cgesvd A double-precision floating-point array for zgesvd	The length of RWORK is LRWORK. Let mx = max(M, N) and mn = min(M, N). If JOBZ = 'N', LRWORK ≥ 5mn. For other values of mx>>mn, LRWORK* ≥ 5mnmn + 5mn. In other cases, LRWORK* ≥ max(5mnmn + 5mn, 2mxmn + 2mn*mn + mn).	Output
iwork	Integer array	Array of dimension 8*min(M,N).	Output
info	Integer	0: The execution is successful. Smaller than 0: If info = -i, the i-th parameter has an illegal value. Greater than 0: An error occurs during computation.	Output

Dependencies

#include "klapack.h"

Examples

C interface:

const char jobz = 'A';
const int m = 5;
const int n = 5;
const int lda = 5;
const int ldu = 5;
const int ldvt = 5;

double a[] = {72.1673, 66.1857, 64.7644, 28.0199, 91.4151,
              6.5180,  62.8483, 72.4323, 46.5760, 8.6928,
              28.9821, 42.1828, 18.6437, 99.8612, 35.6972,
              67.9812, 5.0880,  85.5035, 79.2945, 54.5920,
              28.6869, 49.7512, 7.5186,  28.6929, 84.6041};

double *s = (double*)malloc(m * sizeof(double));
double *u = (double*)malloc(ldu*m * sizeof(double));
double *vt = (double*)malloc(ldvt*n * sizeof(double));
double *iwork = (double*)malloc(8 * m * sizeof(double));

double qwork;
int lwork = -1;
int info = 0;

dgesdd_(&jobz, &m, &n, a, &lda, s, u, &ldu, vt, &ldvt, &qwork, &lwork, iwork, &info);
if (info != 0) {
    printf("Error, info = %d\n", info);
    return info;
}

lwork = (int)qwork;
double *work = (double*)malloc(lwork * sizeof(double));

dgesdd_(&jobz, &m, &n, a, &lda, s, u, &ldu, vt, &ldvt, work, &lwork, iwork, &info);
if (info != 0) {
    printf("Error, info = %d\n", info);
    return info;
}

* output:
* s:
* 254.071516      87.547400       67.835073       64.370735       14.897476
* u:
* -0.395825       -0.377874       -0.458327       -0.488090       -0.502234       -0.155273       -0.217110       0.348128        0.628260        -0.642534       * -0.312281       0.362918        -0.694714       0.528618        0.093312        -0.399962       0.761762        0.321232        -0.289476       -0.269743       * -0.749482       -0.313183       0.287969        0.061841        0.503430
* vt:
* -0.562230       -0.504441       -0.297298       0.148954        -0.564680       -0.340949       0.391046        -0.060652       0.818829        0.238069        * -0.403926       0.372759        0.728616      -0.186514        -0.363629       -0.527964       0.375181        -0.468379       -0.520847       0.299721        * -0.354610       -0.559386       0.397083        -0.035531       0.634351

Fortran interface:

CHARACTER :: jobz = "A"
PARAMETER (m = 5)
PARAMETER (n = 5) 
PARAMETER (lda = 5) 
PARAMETER (ldu = 5)
PARAMETER (ldvt = 5)
INTEGER :: info = 0 
REAL(8) :: a(lad, n)
REAL(8) :: s(m) 
REAL(8) :: u(ldu, m)
REAL(8) :: vt(ldvt, n)
REAL(8) :: iwork(8, m)
REAL(8), ALLOCATABLE :: work(:) 
INTEGER :: lwork = -1 
REAL(8) :: qwork
  
DATA a / 72.1673, 66.1857, 64.7644, 28.0199, 91.4151,
         6.5180,  62.8483, 72.4323, 46.5760, 8.6928,
         28.9821, 42.1828, 18.6437, 99.8612, 35.6972,
         67.9812, 5.0880,  85.5035, 79.2945, 54.5920,
         28.6869, 49.7512, 7.5186,  28.6929, 84.6041/ 
EXTERNAL DGESDD
CALL DGESDD(jobu, jobvt, m, n, a, lda, s, u, ldu, vt, ldvt, qwork, lwork, iwork, info);
lwork = (int)qwork; 
work = (double *)malloc(sizeof(double) * lwork);
CALL DGESDD(jobu, jobvt, m, n, a, lda, s, u, ldu, vt, ldvt, work, lwork, iwork, info);
* 
* Output: 
* s:
* 254.071516      87.547400       67.835073       64.370735       14.897476
* u:
* -0.395825       -0.377874       -0.458327       -0.488090       -0.502234       -0.155273       -0.217110       0.348128        0.628260        -0.642534       * -0.312281       0.362918        -0.694714       0.528618        0.093312        -0.399962       0.761762        0.321232        -0.289476       -0.269743       * -0.749482       -0.313183       0.287969        0.061841        0.503430
* vt:
* -0.562230       -0.504441       -0.297298       0.148954        -0.564680       -0.340949       0.391046        -0.060652       0.818829        0.238069        * -0.403926       0.372759        0.728616      -0.186514        -0.363629       -0.527964       0.375181        -0.468379       -0.520847       0.299721        * -0.354610       -0.559386       0.397083        -0.035531       0.634351

Parent topic: SVD Functions