Transforms the matrix mat
to upper Hessenberg form.
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
integer, | intent(in) | :: | n | Dimension of |
||
integer, | intent(in) | :: | low | First nonzero row |
||
integer, | intent(in) | :: | high | Last nonzero row |
||
real(kind=double), | intent(inout) | :: | mat(0:n-1,0:n-1) | is stored in the lower triangle |
||
integer, | intent(out) | :: | perm(0:n-1) | Permutation vector for elmtrans |
subroutine elmhes(n, low, high, mat, perm )
!+Transforms the matrix `mat` to upper Hessenberg form.
integer,intent(in) :: n
!+Dimension of `mat`
integer,intent(in) :: low
!+First nonzero row
integer,intent(in) :: high
!+Last nonzero row
real(double), intent(inout):: mat(0:n-1,0:n-1)
!+Input: `n`x`n` matrix
!+Output: Upper Hessenberg matrix; additional information on the tranformation
!+is stored in the lower triangle
integer,intent(out) :: perm(0:n-1)
!+Permutation vector for elmtrans
integer :: i, j, m
real(double) :: x, y
do m = low + 1, high-1
i = m
x = ZERO
do j = m, high
if (DABS(mat(j,m-1)) > DABS (x)) then
x = mat(j,m-1)
i = j
endif
enddo
perm(m) = i
if (i.ne.m) then
do j = m - 1, n-1
call RSWAP(mat(i,j), mat(m,j))
enddo
do j = 0, high
call RSWAP(mat(j,i), mat(j,m))
enddo
endif
if (x.ne.ZERO) then
do i = m + 1, high
y = mat(i,m-1)
if (y.ne.ZERO) then
y = y / x
mat(i,m-1) = y
do j = m, n-1
mat(i,j) = mat(i,j) - y * mat(m,j)
enddo
do j = 0, high
mat(j,m) = mat(j,m) + y * mat(j,i)
enddo
endif
enddo !i loop
endif !x <> ZERO
enddo !m loop
end subroutine elmhes