linalg::svd Interface Reference

Computes the singular value decomposition of a matrix A. The SVD is defined as: \( A = U S V^T \), where \( U \) is an M-by-M orthogonal matrix, \( S \) is an M-by-N diagonal matrix, and \( V \) is an N-by-N orthogonal matrix. In the event that \( V \) is complex valued, \( V^H \) is computed instead of \( V^T \). More...

Public Member Functions
	svd_dbl
	svd_cmplx

Detailed Description

Computes the singular value decomposition of a matrix A. The SVD is defined as: \( A = U S V^T \), where \( U \) is an M-by-M orthogonal matrix, \( S \) is an M-by-N diagonal matrix, and \( V \) is an N-by-N orthogonal matrix. In the event that \( V \) is complex valued, \( V^H \) is computed instead of \( V^T \).

Syntax

subroutine svd(real(real64) a(:,:), real(real64) s(:), optional real(real64) u(:,:), optional real(real64) vt(:,:), optional real(real64) work(:), optional integer(int32) olwork, optional class(errors) err)
subroutine svd(complex(real64) a(:,:), real(real64) s(:), optional complex(real64) u(:,:), optional complex(real64) vt(:,:), optional complex(real64) work(:), optional integer(int32) olwork, optional real(real64) rwork(:), optional class(errors) err)

Parameters

[in,out]	a	On input, the M-by-N matrix to factor. The matrix is overwritten on output.
[out]	s	A MIN(M, N)-element array containing the singular values of a sorted in descending order.
[out]	u	An optional argument, that if supplied, is used to contain the orthogonal matrix U from the decomposition. The matrix U contains the left singular vectors, and can be either M-by-M (all left singular vectors are computed), or M-by-MIN(M,N) (only the first MIN(M, N) left singular vectors are computed).
[out]	vt	An optional argument, that if supplied, is used to contain the conjugate transpose of the N-by-N orthogonal matrix V. The matrix V contains the right singular vectors.
[out]	work	An optional input, that if provided, prevents any local memory allocation for complex-valued workspaces. If not provided, the memory required is allocated within. If provided, the length of the array must be at least olwork.
[out]	olwork	An optional output used to determine workspace size. If supplied, the routine determines the optimal size for work, and returns without performing any actual calculations.
[out]	rwork	An optional input, that if provided, prevents any local memory allocation for real-valued workspaces. If not provided, the memory required is allocated within. If provided, the length of the array must be at least 5 * MIN(M, N).
[in,out]	err	An optional errors-based object that if provided can be used to retrieve information relating to any errors encountered during execution. If not provided, a default implementation of the errors class is used internally to provide error handling. Possible errors and warning messages that may be encountered are as follows. LA_ARRAY_SIZE_ERROR: Occurs if any of the input arrays are not sized appropriately. LA_OUT_OF_MEMORY_ERROR: Occurs if local memory must be allocated, and there is insufficient memory available. LA_CONVERGENCE_ERROR: Occurs as a warning if the QR iteration process could not converge to a zero value.

Notes

This routine utilizes the LAPACK routine DGESVD (ZGESVD in the complex case).

See Also

• Wikipedia
• Wolfram MathWorld

Usage

The following example illustrates the calculation of the singular value decomposition of an overdetermined system.

program example
    use iso_fortran_env, only : int32, real64
    use linalg
    implicit none
 
    ! Variables
    real(real64) :: a(3,2), s(2), u(3,3), vt(2,2), ac(3,2)
    integer(int32) :: i
 
    ! Initialize the 3-by-2 matrix A
    !     | 2   1 |
    ! A = |-3   1 |
    !     |-1   1 |
    a = reshape([2.0d0, -3.0d0, -1.0d0, 1.0d0, 1.0d0, 1.0d0], [3, 2])
 
    ! Compute the singular value decomposition of A.  Notice, V**T is returned
    ! instead of V.  Also note, A is overwritten.
    call svd(a, s, u, vt)
 
    ! Display the results
    print '(A)', "U ="
    do i = 1, size(u, 1)
        print *, u(i,:)
    end do
 
    print '(A)', "S ="
    print '(F9.5)', (s(i), i = 1, size(a, 2))
 
    print '(A)', "V**T ="
    do i = 1, size(vt, 1)
        print *, vt(i,:)
    end do
 
    ! Compute U * S * V**T
    call diag_mtx_mult(.true., 1.0d0, s, vt) ! Compute: VT = S * V**T
    ac = matmul(u(:,1:2), vt)
    print '(A)', "U * S * V**T ="
    do i = 1, size(ac, 1)
        print *, ac(i,:)
    end do
end program

The above program produces the following output.

U =
 -0.47411577501825380      -0.81850539032073777      -0.32444284226152509
  0.82566838523833064      -0.28535874325972488      -0.48666426339228758
  0.30575472113569685      -0.49861740208412991       0.81110710565381272
S =
  3.78845
  1.62716
V**T =
 -0.98483334211643059       0.17350299206578967
 -0.17350299206578967      -0.98483334211643059
U * S * V**T =
   1.9999999999999993       0.99999999999999956
  -3.0000000000000000        1.0000000000000000
  -1.0000000000000000       0.99999999999999967

Definition at line 2182 of file linalg.f90.

---

The documentation for this interface was generated from the following file:

• D:/Code/linalg/src/linalg.f90