Balance a real general matrix and isolate eigenvalues whenever possible.
This subroutine is a translation of the ALGOL procedure BALANCE, NUM. MATH. 13, 293-304(1969) by Parlett and Reinsch. HANDBOOK FOR AUTO. COMP., Vol.II-LINEAR ALGEBRA, 315-326(1971).
This subroutine balances a REAL matrix and isolates eigenvalues whenever possible.
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.
A contains the input matrix to be balanced. A is a
two-dimensional REAL array, dimensioned A(NM,N).
A contains the balanced matrix.
LOW and IGH are two INTEGER variables such that A(I,J)
is equal to zero if
(1) I is greater than J and
(2) J=1,...,LOW-1 or I=IGH+1,...,N.
SCALE contains information determining the permutations and
scaling factors used. SCALE is a one-dimensional REAL array,
dimensioned SCALE(N).
Suppose that the principal submatrix in rows LOW through IGH
has been balanced, that P(J) denotes the index interchanged
with J during the permutation step, and that the elements
of the diagonal matrix used are denoted by D(I,J). Then
SCALE(J) = P(J), for J = 1,...,LOW-1
= D(J,J), J = LOW,...,IGH
= P(J) J = IGH+1,...,N.
The order in which the interchanges are made is N to IGH+1,
then 1 TO LOW-1.
Note that 1 is returned for IGH if IGH is zero formally.
The ALGOL procedure EXC contained in BALANCE appears in BALANC in line. (Note that the ALGOL roles of identifiers K,L have been reversed.)
Questions and comments should be directed to B. S. Garbow, Applied Mathematics Division, ARGONNE NATIONAL LABORATORY
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| integer | :: | Nm | ||||
| integer | :: | n | ||||
| real(kind=wp) | :: | a(Nm,*) | ||||
| integer | :: | Low | ||||
| integer | :: | Igh | ||||
| real(kind=wp) | :: | Scale(*) |
subroutine balanc(Nm, n, a, Low, Igh, Scale) implicit none integer :: i, j, k, l, m, n, jj, Nm, Igh, Low, iexc, igo, igo1, igo2 real(wp) :: a(Nm, *), Scale(*) real(wp) :: c, f, g, r, s, b2 logical :: noconv real(wp),parameter :: radix = 16.0_wp b2 = radix*radix k = 1 l = n igo = 1 igo1 = 0 igo2 = 1 ! IN-LINE PROCEDURE FOR ROW AND ! COLUMN EXCHANGE do while (igo2 == 1) igo2 = 0 if (igo1 == 1) then Scale(m) = j if (j /= m) then do i = 1, l f = a(i, j) a(i, j) = a(i, m) a(i, m) = f enddo do i = k, n f = a(j, i) a(j, i) = a(m, i) a(m, i) = f enddo endif if (iexc == 2) then ! SEARCH FOR COLUMNS ISOLATING AN EIGENVALUE ! AND PUSH THEM LEFT k = k + 1 igo = 0 else ! SEARCH FOR ROWS ISOLATING AN EIGENVALUE ! AND PUSH THEM DOWN if (l == 1) then Low = k Igh = l return end if l = l - 1 endif end if ! FOR J=L STEP -1 UNTIL 1 DO -- igo1 = 1 if (igo == 1) then do jj = 1, l igo = 1 j = l + 1 - jj do i = 1, l if (i /= j) then if (a(j, i) /= 0.0_wp) then igo = 0 exit end if endif enddo if (igo == 0) cycle m = l iexc = 1 igo2 = 1 exit enddo if (igo2 == 1) cycle end if do j = k, l igo = 1 do i = k, l if (i /= j) then if (a(i, j) /= 0.0_wp) then igo = 0 exit end if endif enddo if (igo == 0) cycle m = k iexc = 2 igo2 = 1 exit enddo if (igo2 == 1) cycle end do ! NOW BALANCE THE SUBMATRIX IN ROWS K TO L do i = k, l Scale(i) = 1.0_wp enddo ! ITERATIVE LOOP FOR NORM REDUCTION noconv = .true. do while (noconv) noconv = .false. do i = k, l c = 0.0_wp r = 0.0_wp do j = k, l if (j /= i) then c = c + abs(a(j, i)) r = r + abs(a(i, j)) endif enddo ! GUARD AGAINST ZERO C OR R DUE TO UNDERFLOW if (c /= 0.0_wp .and. r /= 0.0_wp) then g = r/radix f = 1.0_wp s = c + r do while (c < g) f = f*radix c = c*b2 end do g = r*radix do while (c >= g) f = f/radix c = c/b2 end do ! NOW BALANCE if ((c + r)/f < 0.95_wp*s) then g = 1.0_wp/f Scale(i) = Scale(i)*f noconv = .true. do j = k, n a(i, j) = a(i, j)*g enddo do j = 1, l a(j, i) = a(j, i)*f enddo endif endif enddo end do Low = k Igh = l end subroutine balanc