gesdd_rel Subroutine

private pure subroutine gesdd_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~~gesdd_rel~~CalledByGraph proc~gesdd_rel forsvd::gesdd_rel proc~svd_rel forsvd::svd_rel proc~svd_rel->proc~gesdd_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 gesdd_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
      integer, dimension(:), allocatable :: iwork     ! Integer work array

      ! External subroutine for calculating the SVD
      interface dgesdd
         pure subroutine dgesdd(f_jobz, f_m, f_n, f_a, f_lda, f_s, f_u, f_ldu, f_vt, f_ldvt, f_work, f_lwork, f_iwork, f_info)
            use kinds
            character,                       intent(in)    :: f_jobz
            integer,                         intent(in)    :: f_m, f_n, f_lda, f_ldu, f_ldvt, f_lwork
            real(rk),  dimension(f_lda, *),  intent(in)    :: f_a
            real(rk),  dimension(*),         intent(out)   :: f_s
            real(rk),  dimension(f_ldu, *),  intent(out)   :: f_u
            real(rk),  dimension(f_ldvt, *), intent(out)   :: f_vt
            real(rk),  dimension(*),         intent(out)   :: f_work
            integer,   dimension(*),         intent(out)   :: f_iwork
            integer,                         intent(out)   :: f_info
         end subroutine dgesdd
      end interface

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

      allocate(iwork(n * 8)) ! Adjust the size as needed

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

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

      deallocate(work)
      deallocate(iwork)
   end subroutine gesdd_rel