svd Interface

public interface svd

Computes the singular value decomposition of a real or complex 2d matrix. (Specification)

Summary

Interface for computing the singular value decomposition of a real or complex 2d matrix.

Description

This interface provides methods for computing the singular value decomposition of a matrix. Supported data types include real and complex. The subroutine returns a real array of singular values, and optionally, left- and right- singular vector matrices, U and V. For a matrix A with size [m,n], full matrix storage for U and V should be [m,m] and [n,n]. It is possible to use partial storage [m,k] and [k,n], k=min(m,n), choosing full_matrices=.false..

Note

The solution is based on LAPACK's singular value decomposition *GESDD methods.

Example

    real(sp) :: a(2,3), s(2), u(2,2), vt(3,3) 
    a = reshape([3,2, 2,3, 2,-2],[2,3])

    call svd(A,s,u,v)
    print *, 'singular values = ',s

Subroutines

private module subroutine stdlib_linalg_svd_c(a, s, u, vt, overwrite_a, full_matrices, err)

Summary

Compute singular value decomposition of a matrix

Description

This function computes the singular value decomposition of a real or complex matrix , and returns the array of singular values, and optionally the left matrix containing the left unitary singular vectors, and the right matrix , containing the right unitary singular vectors.

param: a Input matrix of size [m,n]. param: s Output real array of size [min(m,n)] returning a list of singular values. param: u [optional] Output left singular matrix of size [m,m] or [m,min(m,n)] (.not.full_matrices). Contains singular vectors as columns. param: vt [optional] Output right singular matrix of size [n,n] or [min(m,n),n] (.not.full_matrices). Contains singular vectors as rows.
param: overwrite_a [optional] Flag indicating if the input matrix can be overwritten. param: full_matrices [optional] If .true. (default), matrices and have size [m,m], [n,n]. Otherwise, they are [m,k], [k,n] with k=min(m,n). param: err [optional] State return flag.

Arguments

Type IntentOptional Attributes Name
complex(kind=sp), intent(inout), target :: a(:,:)

Input matrix A[m,n]

real(kind=sp), intent(out) :: s(:)

Array of singular values

complex(kind=sp), intent(out), optional, target :: u(:,:)

The columns of U contain the left singular vectors

complex(kind=sp), intent(out), optional, target :: vt(:,:)

The rows of V^T contain the right singular vectors

logical(kind=lk), intent(in), optional :: overwrite_a

[optional] Can A data be overwritten and destroyed?

logical(kind=lk), intent(in), optional :: full_matrices

[optional] full matrices have shape(u)==[m,m], shape(vh)==[n,n] (default); otherwise they are shape(u)==[m,k] and shape(vh)==[k,n] with k=min(m,n)

type(linalg_state_type), intent(out), optional :: err

[optional] state return flag. On error if not requested, the code will stop

private module subroutine stdlib_linalg_svd_d(a, s, u, vt, overwrite_a, full_matrices, err)

Summary

Compute singular value decomposition of a matrix

Description

This function computes the singular value decomposition of a real or complex matrix , and returns the array of singular values, and optionally the left matrix containing the left unitary singular vectors, and the right matrix , containing the right unitary singular vectors.

param: a Input matrix of size [m,n]. param: s Output real array of size [min(m,n)] returning a list of singular values. param: u [optional] Output left singular matrix of size [m,m] or [m,min(m,n)] (.not.full_matrices). Contains singular vectors as columns. param: vt [optional] Output right singular matrix of size [n,n] or [min(m,n),n] (.not.full_matrices). Contains singular vectors as rows.
param: overwrite_a [optional] Flag indicating if the input matrix can be overwritten. param: full_matrices [optional] If .true. (default), matrices and have size [m,m], [n,n]. Otherwise, they are [m,k], [k,n] with k=min(m,n). param: err [optional] State return flag.

Arguments

Type IntentOptional Attributes Name
real(kind=dp), intent(inout), target :: a(:,:)

Input matrix A[m,n]

real(kind=dp), intent(out) :: s(:)

Array of singular values

real(kind=dp), intent(out), optional, target :: u(:,:)

The columns of U contain the left singular vectors

real(kind=dp), intent(out), optional, target :: vt(:,:)

The rows of V^T contain the right singular vectors

logical(kind=lk), intent(in), optional :: overwrite_a

[optional] Can A data be overwritten and destroyed?

logical(kind=lk), intent(in), optional :: full_matrices

[optional] full matrices have shape(u)==[m,m], shape(vh)==[n,n] (default); otherwise they are shape(u)==[m,k] and shape(vh)==[k,n] with k=min(m,n)

type(linalg_state_type), intent(out), optional :: err

[optional] state return flag. On error if not requested, the code will stop

private module subroutine stdlib_linalg_svd_s(a, s, u, vt, overwrite_a, full_matrices, err)

Summary

Compute singular value decomposition of a matrix

Description

This function computes the singular value decomposition of a real or complex matrix , and returns the array of singular values, and optionally the left matrix containing the left unitary singular vectors, and the right matrix , containing the right unitary singular vectors.

param: a Input matrix of size [m,n]. param: s Output real array of size [min(m,n)] returning a list of singular values. param: u [optional] Output left singular matrix of size [m,m] or [m,min(m,n)] (.not.full_matrices). Contains singular vectors as columns. param: vt [optional] Output right singular matrix of size [n,n] or [min(m,n),n] (.not.full_matrices). Contains singular vectors as rows.
param: overwrite_a [optional] Flag indicating if the input matrix can be overwritten. param: full_matrices [optional] If .true. (default), matrices and have size [m,m], [n,n]. Otherwise, they are [m,k], [k,n] with k=min(m,n). param: err [optional] State return flag.

Arguments

Type IntentOptional Attributes Name
real(kind=sp), intent(inout), target :: a(:,:)

Input matrix A[m,n]

real(kind=sp), intent(out) :: s(:)

Array of singular values

real(kind=sp), intent(out), optional, target :: u(:,:)

The columns of U contain the left singular vectors

real(kind=sp), intent(out), optional, target :: vt(:,:)

The rows of V^T contain the right singular vectors

logical(kind=lk), intent(in), optional :: overwrite_a

[optional] Can A data be overwritten and destroyed?

logical(kind=lk), intent(in), optional :: full_matrices

[optional] full matrices have shape(u)==[m,m], shape(vh)==[n,n] (default); otherwise they are shape(u)==[m,k] and shape(vh)==[k,n] with k=min(m,n)

type(linalg_state_type), intent(out), optional :: err

[optional] state return flag. On error if not requested, the code will stop

private module subroutine stdlib_linalg_svd_z(a, s, u, vt, overwrite_a, full_matrices, err)

Summary

Compute singular value decomposition of a matrix

Description

This function computes the singular value decomposition of a real or complex matrix , and returns the array of singular values, and optionally the left matrix containing the left unitary singular vectors, and the right matrix , containing the right unitary singular vectors.

param: a Input matrix of size [m,n]. param: s Output real array of size [min(m,n)] returning a list of singular values. param: u [optional] Output left singular matrix of size [m,m] or [m,min(m,n)] (.not.full_matrices). Contains singular vectors as columns. param: vt [optional] Output right singular matrix of size [n,n] or [min(m,n),n] (.not.full_matrices). Contains singular vectors as rows.
param: overwrite_a [optional] Flag indicating if the input matrix can be overwritten. param: full_matrices [optional] If .true. (default), matrices and have size [m,m], [n,n]. Otherwise, they are [m,k], [k,n] with k=min(m,n). param: err [optional] State return flag.

Arguments

Type IntentOptional Attributes Name
complex(kind=dp), intent(inout), target :: a(:,:)

Input matrix A[m,n]

real(kind=dp), intent(out) :: s(:)

Array of singular values

complex(kind=dp), intent(out), optional, target :: u(:,:)

The columns of U contain the left singular vectors

complex(kind=dp), intent(out), optional, target :: vt(:,:)

The rows of V^T contain the right singular vectors

logical(kind=lk), intent(in), optional :: overwrite_a

[optional] Can A data be overwritten and destroyed?

logical(kind=lk), intent(in), optional :: full_matrices

[optional] full matrices have shape(u)==[m,m], shape(vh)==[n,n] (default); otherwise they are shape(u)==[m,k] and shape(vh)==[k,n] with k=min(m,n)

type(linalg_state_type), intent(out), optional :: err

[optional] state return flag. On error if not requested, the code will stop