elmhes Subroutine

private subroutine elmhes(Nm, n, Low, Igh, a, Intv)

Reduce a real general matrix to upper Hessenberg form using stabilized elementary similarity transformations.

This subroutine is a translation of the ALGOL procedure ELMHES, NUM. MATH. 12, 349-368(1968) by Martin and Wilkinson. HANDBOOK FOR AUTO. COMP., VOL.II-LINEAR ALGEBRA, 339-358(1971).

Given a REAL GENERAL matrix, this subroutine reduces a submatrix situated in rows and columns LOW through IGH to upper Hessenberg form by stabilized elementary similarity transformations.

On Input

    NM must be set to the row dimension of the two-dimensional
      array parameter, A, as declared in the calling program
      dimension statement.  NM is an INTEGER variable.

    N is the order of the matrix, A.  N is an INTEGER variable.
      N must be less than or equal to NM.

    LOW and IGH are two INTEGER variables determined by the
      balancing subroutine  BALANC.  If  BALANC  has not been
      used, set LOW=1 and IGH equal to the order of the matrix, N.

    A contains the input matrix.  A is a two-dimensional REAL
      array, dimensioned A(NM,N).

On Output

    A contains the upper Hessenberg matrix.  The multipliers which
      were used in the reduction are stored in the remaining
      triangle under the Hessenberg matrix.

    INTV contains information on the rows and columns interchanged
      in the reduction.  Only elements LOW through IGH are used.
      INTV is a one-dimensional INTEGER array, dimensioned INTV(IGH).

 Questions and comments should be directed to B. S. Garbow,
 APPLIED MATHEMATICS DIVISION, ARGONNE NATIONAL LABORATORY

References

  • B. T. Smith, J. M. Boyle, J. J. Dongarra, B. S. Garbow, Y. Ikebe, V. C. Klema and C. B. Moler, Matrix Eigen- system Routines - EISPACK Guide, Springer-Verlag, 1976.

Revision History

  • Author: Smith, B. T., et al.
  • 760101 DATE WRITTEN
  • 890831 Modified array declarations. (WRB)
  • 890831 REVISION DATE from Version 3.2
  • 891214 Prologue converted to Version 4.0 format. (BAB)
  • 920501 Reformatted the REFERENCES section. (WRB)

Arguments

Type IntentOptional Attributes Name
integer :: Nm
integer :: n
integer :: Low
integer :: Igh
real(kind=wp) :: a(Nm,*)
integer :: Intv(*)

Called by

proc~~elmhes~~CalledByGraph proc~elmhes eispack_module::elmhes proc~rg eispack_module::rg proc~rg->proc~elmhes proc~compute_eigenvalues_and_eigenvectors eispack_module::compute_eigenvalues_and_eigenvectors proc~compute_eigenvalues_and_eigenvectors->proc~rg proc~compute_real_eigenvalues_and_normalized_eigenvectors eispack_module::compute_real_eigenvalues_and_normalized_eigenvectors proc~compute_real_eigenvalues_and_normalized_eigenvectors->proc~compute_eigenvalues_and_eigenvectors proc~eispack_test eispack_module::eispack_test proc~eispack_test->proc~compute_eigenvalues_and_eigenvectors

Source Code

 subroutine elmhes(Nm, n, Low, Igh, a, Intv)

    implicit none

    integer :: i, j, m, n, la, Nm, Igh, kp1, Low, mm1, mp1
    real(wp) :: a(Nm, *)
    real(wp) :: x, y
    integer :: Intv(*)

    la = Igh - 1
    kp1 = Low + 1
    if (la >= kp1) then

       do m = kp1, la
          mm1 = m - 1
          x = 0.0_wp
          i = m

          do j = m, Igh
             if (abs(a(j, mm1)) > abs(x)) then
                x = a(j, mm1)
                i = j
             endif
          enddo

          Intv(m) = i
          if (i /= m) then
             ! INTERCHANGE ROWS AND COLUMNS OF A
             do j = mm1, n
                y = a(i, j)
                a(i, j) = a(m, j)
                a(m, j) = y
             enddo

             do j = 1, Igh
                y = a(j, i)
                a(j, i) = a(j, m)
                a(j, m) = y
             enddo
          endif
          ! END INTERCHANGE
          if (x /= 0.0_wp) then
             mp1 = m + 1

             do i = mp1, Igh
                y = a(i, mm1)
                if (y /= 0.0_wp) then
                   y = y/x
                   a(i, mm1) = y

                   do j = m, n
                      a(i, j) = a(i, j) - y*a(m, j)
                   enddo

                   do j = 1, Igh
                      a(j, m) = a(j, m) + y*a(j, i)
                   enddo
                endif

             enddo
          endif

       enddo
    endif

 end subroutine elmhes