1 SUBROUTINE CBAL(NM,N,AR,AI,LOW,IGH,SCALE) 2C***BEGIN PROLOGUE CBAL 3C***DATE WRITTEN 760101 (YYMMDD) 4C***REVISION DATE 830518 (YYMMDD) 5C***CATEGORY NO. D4C1A 6C***KEYWORDS EIGENVALUES,EIGENVECTORS,EISPACK 7C***AUTHOR SMITH, B. T., ET AL. 8C***PURPOSE Balances a complex general matrix and isolates eigenvalues 9C whenever possible. 10C***DESCRIPTION 11C 12C This subroutine is a translation of the ALGOL procedure 13C CBALANCE, which is a complex version of BALANCE, 14C NUM. MATH. 13, 293-304(1969) by Parlett and Reinsch. 15C HANDBOOK FOR AUTO. COMP., VOL.II-LINEAR ALGEBRA, 315-326(1971). 16C 17C This subroutine balances a COMPLEX matrix and isolates 18C eigenvalues whenever possible. 19C 20C On INPUT 21C 22C NM must be set to the row dimension of two-dimensional 23C array parameters as declared in the calling program 24C dimension statement. 25C 26C N is the order of the matrix. 27C 28C AR and AI contain the real and imaginary parts, 29C respectively, of the complex matrix to be balanced. 30C 31C On OUTPUT 32C 33C AR and AI contain the real and imaginary parts, 34C respectively, of the balanced matrix. 35C 36C LOW and IGH are two integers such that AR(I,J) and AI(I,J) 37C are equal to zero if 38C (1) I is greater than J and 39C (2) J=1,...,LOW-1 or I=IGH+1,...,N. 40C 41C SCALE contains information determining the 42C permutations and scaling factors used. 43C 44C Suppose that the principal submatrix in rows LOW through IGH 45C has been balanced, that P(J) denotes the index interchanged 46C with J during the permutation step, and that the elements 47C of the diagonal matrix used are denoted by D(I,J). Then 48C SCALE(J) = P(J), for J = 1,...,LOW-1 49C = D(J,J) J = LOW,...,IGH 50C = P(J) J = IGH+1,...,N. 51C The order in which the interchanges are made is N to IGH+1, 52C then 1 to LOW-1. 53C 54C Note that 1 is returned for IGH if IGH is zero formally. 55C 56C The ALGOL procedure EXC contained in CBALANCE appears in 57C CBAL in line. (Note that the ALGOL roles of identifiers 58C K,L have been reversed.) 59C 60C Questions and comments should be directed to B. S. Garbow, 61C APPLIED MATHEMATICS DIVISION, ARGONNE NATIONAL LABORATORY 62C ------------------------------------------------------------------ 63C***REFERENCES B. T. SMITH, J. M. BOYLE, J. J. DONGARRA, B. S. GARBOW, 64C Y. IKEBE, V. C. KLEMA, C. B. MOLER, *MATRIX EIGEN- 65C SYSTEM ROUTINES - EISPACK GUIDE*, SPRINGER-VERLAG, 66C 1976. 67C***ROUTINES CALLED (NONE) 68C***END PROLOGUE CBAL 69C 70 INTEGER I,J,K,L,M,N,JJ,NM,IGH,LOW,IEXC 71 REAL AR(NM,N),AI(NM,N),SCALE(N) 72 REAL C,F,G,R,S,B2,RADIX 73 LOGICAL NOCONV 74C 75C THE FOLLOWING PORTABLE VALUE OF RADIX WORKS WELL ENOUGH 76C FOR ALL MACHINES WHOSE BASE IS A POWER OF TWO. 77C 78C***FIRST EXECUTABLE STATEMENT CBAL 79 RADIX = 16 80C 81 B2 = RADIX * RADIX 82 K = 1 83 L = N 84 GO TO 100 85C .......... IN-LINE PROCEDURE FOR ROW AND 86C COLUMN EXCHANGE .......... 87 20 SCALE(M) = J 88 IF (J .EQ. M) GO TO 50 89C 90 DO 30 I = 1, L 91 F = AR(I,J) 92 AR(I,J) = AR(I,M) 93 AR(I,M) = F 94 F = AI(I,J) 95 AI(I,J) = AI(I,M) 96 AI(I,M) = F 97 30 CONTINUE 98C 99 DO 40 I = K, N 100 F = AR(J,I) 101 AR(J,I) = AR(M,I) 102 AR(M,I) = F 103 F = AI(J,I) 104 AI(J,I) = AI(M,I) 105 AI(M,I) = F 106 40 CONTINUE 107C 108 50 GO TO (80,130), IEXC 109C .......... SEARCH FOR ROWS ISOLATING AN EIGENVALUE 110C AND PUSH THEM DOWN .......... 111 80 IF (L .EQ. 1) GO TO 280 112 L = L - 1 113C .......... FOR J=L STEP -1 UNTIL 1 DO -- .......... 114 100 DO 120 JJ = 1, L 115 J = L + 1 - JJ 116C 117 DO 110 I = 1, L 118 IF (I .EQ. J) GO TO 110 119 IF (AR(J,I) .NE. 0.0E0 .OR. AI(J,I) .NE. 0.0E0) GO TO 120 120 110 CONTINUE 121C 122 M = L 123 IEXC = 1 124 GO TO 20 125 120 CONTINUE 126C 127 GO TO 140 128C .......... SEARCH FOR COLUMNS ISOLATING AN EIGENVALUE 129C AND PUSH THEM LEFT .......... 130 130 K = K + 1 131C 132 140 DO 170 J = K, L 133C 134 DO 150 I = K, L 135 IF (I .EQ. J) GO TO 150 136 IF (AR(I,J) .NE. 0.0E0 .OR. AI(I,J) .NE. 0.0E0) GO TO 170 137 150 CONTINUE 138C 139 M = K 140 IEXC = 2 141 GO TO 20 142 170 CONTINUE 143C .......... NOW BALANCE THE SUBMATRIX IN ROWS K TO L .......... 144 DO 180 I = K, L 145 180 SCALE(I) = 1.0E0 146C .......... ITERATIVE LOOP FOR NORM REDUCTION .......... 147 190 NOCONV = .FALSE. 148C 149 DO 270 I = K, L 150 C = 0.0E0 151 R = 0.0E0 152C 153 DO 200 J = K, L 154 IF (J .EQ. I) GO TO 200 155 C = C + ABS(AR(J,I)) + ABS(AI(J,I)) 156 R = R + ABS(AR(I,J)) + ABS(AI(I,J)) 157 200 CONTINUE 158C .......... GUARD AGAINST ZERO C OR R DUE TO UNDERFLOW .......... 159 IF (C .EQ. 0.0E0 .OR. R .EQ. 0.0E0) GO TO 270 160 G = R / RADIX 161 F = 1.0E0 162 S = C + R 163 210 IF (C .GE. G) GO TO 220 164 F = F * RADIX 165 C = C * B2 166 GO TO 210 167 220 G = R * RADIX 168 230 IF (C .LT. G) GO TO 240 169 F = F / RADIX 170 C = C / B2 171 GO TO 230 172C .......... NOW BALANCE .......... 173 240 IF ((C + R) / F .GE. 0.95E0 * S) GO TO 270 174 G = 1.0E0 / F 175 SCALE(I) = SCALE(I) * F 176 NOCONV = .TRUE. 177C 178 DO 250 J = K, N 179 AR(I,J) = AR(I,J) * G 180 AI(I,J) = AI(I,J) * G 181 250 CONTINUE 182C 183 DO 260 J = 1, L 184 AR(J,I) = AR(J,I) * F 185 AI(J,I) = AI(J,I) * F 186 260 CONTINUE 187C 188 270 CONTINUE 189C 190 IF (NOCONV) GO TO 190 191C 192 280 LOW = K 193 IGH = L 194 RETURN 195 END 196