stdlib_linalg.fypp Source File


This file depends on

sourcefile~~stdlib_linalg.fypp~~EfferentGraph sourcefile~stdlib_linalg.fypp stdlib_linalg.fypp sourcefile~stdlib_kinds.f90 stdlib_kinds.f90 sourcefile~stdlib_linalg.fypp->sourcefile~stdlib_kinds.f90

Files dependent on this one

sourcefile~~stdlib_linalg.fypp~~AfferentGraph sourcefile~stdlib_linalg.fypp stdlib_linalg.fypp sourcefile~stdlib_stats_corr.fypp stdlib_stats_corr.fypp sourcefile~stdlib_stats_corr.fypp->sourcefile~stdlib_linalg.fypp sourcefile~stdlib_linalg_diag.fypp stdlib_linalg_diag.fypp sourcefile~stdlib_linalg_diag.fypp->sourcefile~stdlib_linalg.fypp

Contents

Source Code


Source Code

#:include "common.fypp"
#:set RCI_KINDS_TYPES = REAL_KINDS_TYPES + CMPLX_KINDS_TYPES + INT_KINDS_TYPES
module stdlib_linalg
  !!Provides a support for various linear algebra procedures
  !! ([Specification](../page/specs/stdlib_linalg.html))
  use stdlib_kinds, only: sp, dp, qp, &
    int8, int16, int32, int64
  implicit none
  private

  public :: diag
  public :: eye
  public :: trace

  interface diag
    !! version: experimental
    !!
    !! Creates a diagonal array or extract the diagonal elements of an array
    !! ([Specification](../page/specs/stdlib_linalg.html#description))
      !
      ! Vector to matrix
      !
    #:for k1, t1 in RCI_KINDS_TYPES
      module function diag_${t1[0]}$${k1}$(v) result(res)
        ${t1}$, intent(in) :: v(:)
        ${t1}$ :: res(size(v),size(v))
      end function diag_${t1[0]}$${k1}$
    #:endfor
    #:for k1, t1 in RCI_KINDS_TYPES
      module function diag_${t1[0]}$${k1}$_k(v,k) result(res)
        ${t1}$, intent(in) :: v(:)
        integer, intent(in) :: k
        ${t1}$ :: res(size(v)+abs(k),size(v)+abs(k))
      end function diag_${t1[0]}$${k1}$_k
    #:endfor

      !
      ! Matrix to vector
      !
    #:for k1, t1 in RCI_KINDS_TYPES
      module function diag_${t1[0]}$${k1}$_mat(A) result(res)
        ${t1}$, intent(in) :: A(:,:)
        ${t1}$ :: res(minval(shape(A)))
      end function diag_${t1[0]}$${k1}$_mat
    #:endfor
    #:for k1, t1 in RCI_KINDS_TYPES
      module function diag_${t1[0]}$${k1}$_mat_k(A,k) result(res)
        ${t1}$, intent(in) :: A(:,:)
        integer, intent(in) :: k
        ${t1}$ :: res(minval(shape(A))-abs(k))
      end function diag_${t1[0]}$${k1}$_mat_k
    #:endfor
  end interface

  ! Matrix trace
  interface trace
    !! version: experimental
    !!
    !! Computes the trace of a matrix
    !! ([Specification](../page/specs/stdlib_linalg.html#description_2))
    #:for k1, t1 in RCI_KINDS_TYPES
      module procedure trace_${t1[0]}$${k1}$
    #:endfor
  end interface

contains

    function eye(n) result(res)
      !! version: experimental
      !!
      !! Constructs the identity matrix
      !! ([Specification](../page/specs/stdlib_linalg.html#description_1))
      integer, intent(in) :: n
      integer(int8) :: res(n, n)
      integer :: i
      res = 0
      do i = 1, n
         res(i, i) = 1
      end do
    end function eye


    #:for k1, t1 in RCI_KINDS_TYPES
      function trace_${t1[0]}$${k1}$(A) result(res)
        ${t1}$, intent(in) :: A(:,:)
        ${t1}$ :: res
        integer :: i
        res = 0
        do i = 1, minval(shape(A))
          res = res + A(i,i)
        end do
      end function trace_${t1[0]}$${k1}$
    #:endfor
end module