gesvd_rel Subroutine

private pure subroutine gesvd_rel(A, U, S, VT)

Calculates the singular value decomposition (SVD) of a matrix A.

Arguments

Type IntentOptional Attributes Name
real(kind=rk), intent(in), dimension(:, :), contiguous :: A
real(kind=rk), intent(out), dimension(:,:) :: U
real(kind=rk), intent(out), dimension(:) :: S
real(kind=rk), intent(out), dimension(:,:) :: VT

Called by

proc~~gesvd_rel~~CalledByGraph proc~gesvd_rel forsvd::gesvd_rel proc~svd_rel forsvd::svd_rel proc~svd_rel->proc~gesvd_rel interface~svd forsvd::svd interface~svd->proc~svd_rel proc~lowrank forsvd::tsvd%lowrank proc~lowrank->interface~svd program~benchmark benchmark program~benchmark->interface~svd program~test1 test1 program~test1->interface~svd program~test2 test2 program~test2->interface~svd program~test3 test3 program~test3->interface~svd proc~compress_pixels forsvd::pixel%compress_pixels proc~compress_pixels->proc~lowrank program~test4 test4 program~test4->proc~lowrank program~test5 test5 program~test5->proc~compress_pixels

Source Code

   pure subroutine gesvd_rel(A, U,S,VT)

      ! Inputs:
      real(rk), dimension(:, :), contiguous,          intent(in)  :: A    ! Input matrix A

      ! Outputs:
      real(rk), dimension(:,:), intent(out) :: U    ! Left singular vectors
      real(rk), dimension(:,:), intent(out) :: VT   ! Right singular vectors
      real(rk), dimension(:),   intent(out) :: S    ! Singular values

      ! Local variables
      real(rk)                                                    :: work1(1) ! memory allocation query
      real(rk), dimension(:), allocatable                         :: work     ! Work array
      integer                                                     :: m, n, lwork, info

      ! External subroutine for calculating the SVD
      interface dgesvd
         pure subroutine dgesvd(jobuf,jobvtf,mf,nf,af,ldaf,sf,uf,lduf,vtf,ldvtf,workf,lworkf,infof)
            use kinds
            character, intent(in)  :: jobuf, jobvtf
            integer,   intent(in)  :: mf, nf, ldaf, lduf, ldvtf, lworkf
            real(rk),  intent(in)  :: Af(ldaf, *)
            real(rk),  intent(out) :: Sf(min(mf, nf))
            real(rk),  intent(out) :: Uf(lduf, *), VTf(ldvtf, *)
            real(rk),  intent(out) :: workf(*)
            integer,   intent(out) :: infof
         end subroutine dgesvd
      end interface

      m = size(A, 1)
      n = size(A, 2)

      ! Calculate the optimal size of the work array
      call dgesvd('S', 'S', m, n, A, m, S, U, m, VT, n, work1, -1, info)
      lwork = nint(work1(1))
      allocate(work(lwork))

      call dgesvd('S', 'S', m, n, A, m, S, U, m, VT, n, work, lwork, info)

      deallocate(work)
   end subroutine gesvd_rel