invert_matrix Function

public function invert_matrix(A, errCode) result(Ainv)

Returns the inverse of a matrix calculated by finding the LU decomposition. Depends on LAPACK.

Arguments

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

Matrix to invert
integer, intent(out), optional :: errCode

If error code is present return error code and do not abort

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

inverse of A


Calls

proc~~invert_matrix~~CallsGraph proc~invert_matrix invert_matrix dgetrf dgetrf proc~invert_matrix->dgetrf dgetri dgetri proc~invert_matrix->dgetri proc~tem_abort tem_abort proc~invert_matrix->proc~tem_abort mpi_abort mpi_abort proc~tem_abort->mpi_abort

Called by

proc~~invert_matrix~~CalledByGraph proc~invert_matrix invert_matrix proc~build_matrixlsf_linearintp build_matrixLSF_linearIntp proc~build_matrixlsf_linearintp->proc~invert_matrix proc~build_matrixlsf_quadintp build_matrixLSF_quadIntp proc~build_matrixlsf_quadintp->proc~invert_matrix proc~append_intpmatrixlsf append_intpMatrixLSF proc~append_intpmatrixlsf->proc~build_matrixlsf_linearintp proc~append_intpmatrixlsf->proc~build_matrixlsf_quadintp interface~append~9 append interface~append~9->proc~append_intpmatrixlsf

Source Code

  function invert_matrix(A, errCode) result(Ainv)
    ! ---------------------------------------------------------------------------
    !> Matrix to invert
    real(kind=rk), dimension(:,:), intent(in) :: A
    !> If error code is present return error code and do not abort
    integer, optional, intent(out) :: errCode
    !> inverse of A
    real(kind=rk), dimension(size(A,1),size(A,2)) :: Ainv
    ! ---------------------------------------------------------------------------
    real(kind=rk), dimension(size(A,1)) :: work  ! work array for LAPACK
    integer, dimension(size(A,1)) :: ipiv   ! pivot indices
    integer :: n, info
    ! External procedures defined in LAPACK
    external DGETRF
    external DGETRI
    ! ---------------------------------------------------------------------------
    ! Store A in Ainv to prevent it from being overwritten by LAPACK
    Ainv = A
    n = size(A,1)

    ! DGETRF computes an LU factorization of a general M-by-N matrix A
    ! using partial pivoting with row interchanges.
    call DGETRF(n, n, Ainv, n, ipiv, info)

    if (info /= 0) then
      if (present(errCode)) then
        errCode = info
        write(logUnit(5),*) 'WARNING: Matrix is numerically singular!'
        return
      else
        call tem_abort('Matrix is numerically singular!')
      end if
    end if

    ! DGETRI computes the inverse of a matrix using the LU factorization
    ! computed by DGETRF.
    call DGETRI(n, Ainv, n, ipiv, work, n, info)

    if (info /= 0) then
      if (present(errCode)) then
        errCode = info
        write(logUnit(5),*) 'WARNING: Matrix inversion failed!'
        return
      else
        call tem_abort('Matrix inversion failed!')
      end if
    end if

    ! successfull
    if (present(errCode)) errCode = 0
  end function invert_matrix