SVD computation?

[Christian Waechter]

Hi Blaze-team,

I went through the blaze-wiki and could not find any description to a SVD routine in the lapack section. Do you offer such a computation and I just could not find it?

In case you do not have it yet, are there plans to integrate it into the library?

Cheers, Christian

[Klaus Iglberger]

Hi Christian!

Thanks a lot for the proposal. Unfortunately SVD is not yet available, but we definitely plan to integrate it. We will keep the issue open until the feature is available.

Best regards,

Klaus!

[Klaus Iglberger]

Summary

The feature has been implemented, tested, optimized, and documented as required. It is immediately available via cloning the Blaze repository and will be officially released in Blaze 3.1.

Singular Value Decomposition

The singular value decomposition (SVD) of a dense matrix can be computed via the svd() functions:

template< typename MT, bool SO, typename VT, bool TF >
void svd( const DenseMatrix<MT,SO>& A, DenseVector<VT,TF>& s );

template< typename MT1, bool SO, typename VT, bool TF, typename MT2, typename MT3 >
void svd( const DenseMatrix<MT1,SO>& A, DenseMatrix<MT2,SO>& U, DenseVector<VT,TF>& s, DenseMatrix<MT3,SO>& V );

template< typename MT, bool SO, typename VT, bool TF, typename ST >
size_t svd( const DenseMatrix<MT,SO>& A, DenseVector<VT,TF>& s, ST low, ST upp );

template< typename MT1, bool SO, typename VT, bool TF, typename MT2, typename MT3, typename ST >
size_t svd( const DenseMatrix<MT1,SO>& A, DenseMatrix<MT2,SO>& U, DenseVector<VT,TF>& s, DenseMatrix<MT3,SO>& V, ST low, ST upp );

The first and third function compute only singular values of the given general m-by-n matrix, the second and fourth function additionally compute singular vectors. The resulting singular values are returned in the given vector s, the left singular vectors are returned in the given matrix U, and the right singular vectors are returned in the matrix V. s, U, and V are resized to the correct dimensions (if possible and necessary).

The third and fourth function allow for the specification of a subset of singular values and/or vectors. The number of singular values and vectors to be computed is specified by the lower bound low and the upper bound upp, which either form an integral or a floating point range.

In case low and upp form are of integral type, the function computes all singular values in the index range [low..upp]. The num resulting real and non-negative singular values are stored in descending order in the given vector s, which is either resized (if possible) or expected to be a num-dimensional vector. The resulting left singular vectors are stored in the given matrix U, which is either resized (if possible) or expected to be a m-by-num matrix. The resulting right singular vectors are stored in the given matrix V, which is either resized (if possible) or expected to be a num-by-n matrix.

In case low and upp are of floating point type, the function computes all singular values in the half-open interval (low..upp]. The resulting real and non-negative singular values are stored in descending order in the given vector s, which is either resized (if possible) or expected to be a min(m,n)-dimensional vector. The resulting left singular vectors are stored in the given matrix U, which is either resized (if possible) or expected to be a m-by-min(m,n) matrix. The resulting right singular vectors are stored in the given matrix V, which is either resized (if possible) or expected to be a min(m,n)-by-n matrix.

The functions fail if ...

• ... the given matrix U is a fixed size matrix and the dimensions don't match;
• ... the given vector s is a fixed size vector and the size doesn't match;
• ... the given matrix V is a fixed size matrix and the dimensions don't match;
• ... the given scalar values don't form a proper range;
• ... the singular value decomposition fails.

In all failure cases an exception is thrown.

Examples:

using blaze::DynamicMatrix;
using blaze::DynamicVector;
using blaze::rowMajor;
using blaze::columnVector;

DynamicMatrix<double,rowMajor>  A( 5UL, 8UL );  // The general matrix A
// ... Initialization

DynamicMatrix<double,rowMajor>     U;  // The matrix for the left singular vectors
DynamicVector<double,columnVector> s;  // The vector for the singular values
DynamicMatrix<double,rowMajor>     V;  // The matrix for the right singular vectors

svd( A, U, s, V );

using blaze::DynamicMatrix;
using blaze::DynamicVector;
using blaze::rowMajor;
using blaze::columnVector;

DynamicMatrix<complex<double>,rowMajor>  A( 5UL, 8UL );  // The general matrix A
// ... Initialization

DynamicMatrix<complex<double>,rowMajor> U;  // The matrix for the left singular vectors
DynamicVector<double,columnVector>      s;  // The vector for the singular values
DynamicMatrix<complex<double>,rowMajor> V;  // The matrix for the right singular vectors

svd( A, U, s, V, 0, 2 );

Note that all svd() functions can only be used for dense matrices with float, double, complex<float> or complex<double> element type. The attempt to call the function with matrices of any other element type or with a sparse matrix results in a compile time error!

Also note that the functions compute the singular values and/or singular vectors of a dense matrix by means of LAPACK kernels. Thus the functions can only be used if the fitting LAPACK library is available and linked to the executable. Otherwise a linker error will be created.