Module matrix_mod

module matrix_mod

        ! Uses
    use precision_mod

        ! Variables
    logical, public :: MatrixFail

        ! Interfaces
    public interface OPERATOR (.TIMES.)
    public interface OPERATOR (.TRANSPOSE.)
    public interface OPERATOR (.INVERT.)

        ! Subroutines and functions
    public function R4MatrixTranspose (A) result (Atranspose)
    public function R8MatrixTranspose (A) result (Atranspose)
    public function R4MatrixMultiply (A, B) result (AB)
    public function R8MatrixMultiply (A, B) result (AB)
    public function R4MatrixInvert (A) result (Ainverse)
    public function R8MatrixInvert (A) result (Ainverse)
    public function MatSolv (A, B, x) result (iStat)
    public recursive subroutine Triangulate (Aug, iStat)

end module matrix_mod

Description of Variables

logical, public :: MatrixFail

Description of Interfaces

public interface OPERATOR (.TIMES.)
    module procedure R4MatrixMultiply
    module procedure R8MatrixMultiply
end interface OPERATOR (.TIMES.)

public interface OPERATOR (.TRANSPOSE.)
    module procedure R4MatrixTranspose
    module procedure R8MatrixTranspose
end interface OPERATOR (.TRANSPOSE.)

public interface OPERATOR (.INVERT.)
    module procedure R4MatrixInvert
    module procedure R8MatrixInvert
end interface OPERATOR (.INVERT.)

Description of Subroutines and Functions

public function R4MatrixTranspose (A) result (Atranspose)
    real (kind=R4), DIMENSION(:,:), INTENT(IN) :: A
    real (kind=R4), DIMENSION(SIZE(A,2),SIZE(A,1)) :: Atranspose
end function R4MatrixTranspose

public function R8MatrixTranspose (A) result (Atranspose)
    real (kind=R8), DIMENSION(:,:), INTENT(IN) :: A
    real (kind=R8), DIMENSION(SIZE(A,2),SIZE(A,1)) :: Atranspose
end function R8MatrixTranspose

public function R4MatrixMultiply (A, B) result (AB)
    real (kind=R4), DIMENSION(:,:), INTENT(IN) :: A
    real (kind=R4), DIMENSION(:,:), INTENT(IN) :: B
    real (kind=R4), DIMENSION(SIZE(A,1),SIZE(B,2)) :: AB
end function R4MatrixMultiply

public function R8MatrixMultiply (A, B) result (AB)
    real (kind=R8), DIMENSION(:,:), INTENT(IN) :: A
    real (kind=R8), DIMENSION(:,:), INTENT(IN) :: B
    real (kind=R8), DIMENSION(SIZE(A,1),SIZE(B,2)) :: AB
end function R8MatrixMultiply

public function R4MatrixInvert (A) result (Ainverse)
    real (kind=R4), DIMENSION(:,:), INTENT(IN) :: A
    real (kind=R4), DIMENSION(SIZE(A,1),SIZE(A,2)) :: Ainverse
    ! Calls: kerror
end function R4MatrixInvert

public function R8MatrixInvert (A) result (Ainverse)
    real (kind=R8), DIMENSION(:,:), INTENT(IN) :: A
    real (kind=R8), DIMENSION(SIZE(A,1),SIZE(A,2)) :: Ainverse
    ! Calls: kerror
end function R8MatrixInvert

public function MatSolv (A, B, x) result (iStat)
    real (kind=R8), DIMENSION(:, :), INTENT(in) :: A
    real (kind=R8), DIMENSION(:), INTENT(in) :: B
    real (kind=R8), DIMENSION(SIZE(A, 2)), INTENT(out) :: x
    integer (kind=I4) :: iStat
    ! Calls: Triangulate
end function MatSolv
$  SUBROUTINE test_matrix
$
$    REAL, DIMENSION(3,3):: B
$    REAL, DIMENSION(SIZE(B,2),SIZE(B,1)):: Btranspose
$    REAL, DIMENSION(SIZE(B,1),SIZE(B,2)):: Binverse
$    REAL, DIMENSION(SIZE(B,1),SIZE(B,2)):: BBinv
$    
$    B(1,1) = 1; B(1,2) = 2; B(1,3) = 5
$    B(2,1) = 2; B(2,2) = 3; B(2,3) = 8
$    B(3,1) =-1; B(3,2) = 1; B(3,3) = 2
$    
$    Btranspose = .transpose.B
$    
$    print*, 'Transpose:'
$    Print*, Btranspose(1,1), Btranspose(1,2), Btranspose(1,3)
$    Print*, Btranspose(2,1), Btranspose(2,2), Btranspose(2,3)
$    Print*, Btranspose(3,1), Btranspose(3,2), Btranspose(3,3)
$    
$    Binverse = .invert.B
$    
$    print*, 'Inverse:'
$    Print*, Binverse(1,1), Binverse(1,2), Binverse(1,3)
$    Print*, Binverse(2,1), Binverse(2,2), Binverse(2,3)
$    Print*, Binverse(3,1), Binverse(3,2), Binverse(3,3)
$    
$    BBinv = B.times.(.invert.B)
$    print*, 'Mult'
$    Print*, BBinv(1,1), BBinv(1,2), BBinv(1,3)
$    Print*, BBinv(2,1), BBinv(2,2), BBinv(2,3)
$    Print*, BBinv(3,1), BBinv(3,2), BBinv(3,3)
$    
$  END SUBROUTINE test_matrix
------------------------------------------------------------------------------
 MatSolv solves a matrix equation of the type Ax = B.

public recursive subroutine Triangulate (Aug, iStat)
    real (kind=R8), DIMENSION(:, :), INTENT(inout) :: Aug
    integer (kind=I4), INTENT(inout) :: iStat
    ! Calls: Triangulate
end subroutine Triangulate
------------------------------------------------------------------------------
 Triangulate produces a triangular augmented matrix equivalent to the input
 augmented matrix.
------------------------------------------------------------------------------