gemm Function

private pure function gemm(A, B) result(C)

Calculates the matrix-matrix product of two matrices A and B.

Arguments

Type IntentOptional Attributes Name
real(kind=rk), intent(in), dimension(:, :), contiguous :: A
real(kind=rk), intent(in), dimension(:, :), contiguous :: B

Return Value real(kind=rk), dimension(size(A,1), size(B,2))


Called by

proc~~gemm~~CalledByGraph proc~gemm forinv::gemm proc~pinvlu_rel forinv::pinvLU_rel proc~pinvlu_rel->proc~gemm proc~pinv_rel forinv::pinv_rel proc~pinv_rel->proc~pinvlu_rel interface~inv forinv::inv interface~inv->proc~pinv_rel program~test1 test1 program~test1->interface~inv program~test2 test2 program~test2->interface~inv program~test3 test3 program~test3->interface~inv

Source Code

   pure function gemm(A, B) result(C)
#elif defined (IMPURE)
   impure function gemm(A, B) result(C)
#endif

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

      ! Outputs:
      real(rk), dimension(size(A,1), size(B,2))          :: C     ! Matrix-matrix product of A and B

      ! Local variables
      integer                                            :: m, n, k

      ! External subroutine for calculating the matrix-matrix product of two matrices A and B.
      interface
#if defined (PURE)
         pure subroutine dgemm(f_transa, f_transb, f_m, f_n, f_k, f_alpha, f_a, f_lda, f_b, f_ldb, f_beta, f_c, f_ldc)
#elif defined (IMPURE)
         impure subroutine dgemm(f_transa, f_transb, f_m, f_n, f_k, f_alpha, f_a, f_lda, f_b, f_ldb, f_beta, f_c, f_ldc)
#endif
            import rk
            integer,   intent(in)    :: f_ldc
            integer,   intent(in)    :: f_ldb
            integer,   intent(in)    :: f_lda
            character, intent(in)    :: f_transa
            character, intent(in)    :: f_transb
            integer,   intent(in)    :: f_m
            integer,   intent(in)    :: f_n
            integer,   intent(in)    :: f_k
            real(rk),  intent(in)    :: f_alpha
            real(rk),  intent(in)    :: f_a(f_lda, *)
            real(rk),  intent(in)    :: f_b(f_ldb, *)
            real(rk),  intent(in)    :: f_beta
            real(rk),  intent(inout) :: f_c(f_ldc, *)
         end subroutine dgemm
      end interface

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

      call dgemm('N', 'N', m, n, k, 1.0_rk, A, m, B, k, 0.0_rk, C, m)

   end function gemm