001:       SUBROUTINE CGBEQU( M, N, KL, KU, AB, LDAB, R, C, ROWCND, COLCND,
002:      $                   AMAX, INFO )
003: *
004: *  -- LAPACK routine (version 3.2) --
005: *     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
006: *     November 2006
007: *
008: *     .. Scalar Arguments ..
009:       INTEGER            INFO, KL, KU, LDAB, M, N
010:       REAL               AMAX, COLCND, ROWCND
011: *     ..
012: *     .. Array Arguments ..
013:       REAL               C( * ), R( * )
014:       COMPLEX            AB( LDAB, * )
015: *     ..
016: *
017: *  Purpose
018: *  =======
019: *
020: *  CGBEQU computes row and column scalings intended to equilibrate an
021: *  M-by-N band matrix A and reduce its condition number.  R returns the
022: *  row scale factors and C the column scale factors, chosen to try to
023: *  make the largest element in each row and column of the matrix B with
024: *  elements B(i,j)=R(i)*A(i,j)*C(j) have absolute value 1.
025: *
026: *  R(i) and C(j) are restricted to be between SMLNUM = smallest safe
027: *  number and BIGNUM = largest safe number.  Use of these scaling
028: *  factors is not guaranteed to reduce the condition number of A but
029: *  works well in practice.
030: *
031: *  Arguments
032: *  =========
033: *
034: *  M       (input) INTEGER
035: *          The number of rows of the matrix A.  M >= 0.
036: *
037: *  N       (input) INTEGER
038: *          The number of columns of the matrix A.  N >= 0.
039: *
040: *  KL      (input) INTEGER
041: *          The number of subdiagonals within the band of A.  KL >= 0.
042: *
043: *  KU      (input) INTEGER
044: *          The number of superdiagonals within the band of A.  KU >= 0.
045: *
046: *  AB      (input) COMPLEX array, dimension (LDAB,N)
047: *          The band matrix A, stored in rows 1 to KL+KU+1.  The j-th
048: *          column of A is stored in the j-th column of the array AB as
049: *          follows:
050: *          AB(ku+1+i-j,j) = A(i,j) for max(1,j-ku)<=i<=min(m,j+kl).
051: *
052: *  LDAB    (input) INTEGER
053: *          The leading dimension of the array AB.  LDAB >= KL+KU+1.
054: *
055: *  R       (output) REAL array, dimension (M)
056: *          If INFO = 0, or INFO > M, R contains the row scale factors
057: *          for A.
058: *
059: *  C       (output) REAL array, dimension (N)
060: *          If INFO = 0, C contains the column scale factors for A.
061: *
062: *  ROWCND  (output) REAL
063: *          If INFO = 0 or INFO > M, ROWCND contains the ratio of the
064: *          smallest R(i) to the largest R(i).  If ROWCND >= 0.1 and
065: *          AMAX is neither too large nor too small, it is not worth
066: *          scaling by R.
067: *
068: *  COLCND  (output) REAL
069: *          If INFO = 0, COLCND contains the ratio of the smallest
070: *          C(i) to the largest C(i).  If COLCND >= 0.1, it is not
071: *          worth scaling by C.
072: *
073: *  AMAX    (output) REAL
074: *          Absolute value of largest matrix element.  If AMAX is very
075: *          close to overflow or very close to underflow, the matrix
076: *          should be scaled.
077: *
078: *  INFO    (output) INTEGER
079: *          = 0:  successful exit
080: *          < 0:  if INFO = -i, the i-th argument had an illegal value
081: *          > 0:  if INFO = i, and i is
082: *                <= M:  the i-th row of A is exactly zero
083: *                >  M:  the (i-M)-th column of A is exactly zero
084: *
085: *  =====================================================================
086: *
087: *     .. Parameters ..
088:       REAL               ONE, ZERO
089:       PARAMETER          ( ONE = 1.0E+0, ZERO = 0.0E+0 )
090: *     ..
091: *     .. Local Scalars ..
092:       INTEGER            I, J, KD
093:       REAL               BIGNUM, RCMAX, RCMIN, SMLNUM
094:       COMPLEX            ZDUM
095: *     ..
096: *     .. External Functions ..
097:       REAL               SLAMCH
098:       EXTERNAL           SLAMCH
099: *     ..
100: *     .. External Subroutines ..
101:       EXTERNAL           XERBLA
102: *     ..
103: *     .. Intrinsic Functions ..
104:       INTRINSIC          ABS, AIMAG, MAX, MIN, REAL
105: *     ..
106: *     .. Statement Functions ..
107:       REAL               CABS1
108: *     ..
109: *     .. Statement Function definitions ..
110:       CABS1( ZDUM ) = ABS( REAL( ZDUM ) ) + ABS( AIMAG( ZDUM ) )
111: *     ..
112: *     .. Executable Statements ..
113: *
114: *     Test the input parameters
115: *
116:       INFO = 0
117:       IF( M.LT.0 ) THEN
118:          INFO = -1
119:       ELSE IF( N.LT.0 ) THEN
120:          INFO = -2
121:       ELSE IF( KL.LT.0 ) THEN
122:          INFO = -3
123:       ELSE IF( KU.LT.0 ) THEN
124:          INFO = -4
125:       ELSE IF( LDAB.LT.KL+KU+1 ) THEN
126:          INFO = -6
127:       END IF
128:       IF( INFO.NE.0 ) THEN
129:          CALL XERBLA( 'CGBEQU', -INFO )
130:          RETURN
131:       END IF
132: *
133: *     Quick return if possible
134: *
135:       IF( M.EQ.0 .OR. N.EQ.0 ) THEN
136:          ROWCND = ONE
137:          COLCND = ONE
138:          AMAX = ZERO
139:          RETURN
140:       END IF
141: *
142: *     Get machine constants.
143: *
144:       SMLNUM = SLAMCH( 'S' )
145:       BIGNUM = ONE / SMLNUM
146: *
147: *     Compute row scale factors.
148: *
149:       DO 10 I = 1, M
150:          R( I ) = ZERO
151:    10 CONTINUE
152: *
153: *     Find the maximum element in each row.
154: *
155:       KD = KU + 1
156:       DO 30 J = 1, N
157:          DO 20 I = MAX( J-KU, 1 ), MIN( J+KL, M )
158:             R( I ) = MAX( R( I ), CABS1( AB( KD+I-J, J ) ) )
159:    20    CONTINUE
160:    30 CONTINUE
161: *
162: *     Find the maximum and minimum scale factors.
163: *
164:       RCMIN = BIGNUM
165:       RCMAX = ZERO
166:       DO 40 I = 1, M
167:          RCMAX = MAX( RCMAX, R( I ) )
168:          RCMIN = MIN( RCMIN, R( I ) )
169:    40 CONTINUE
170:       AMAX = RCMAX
171: *
172:       IF( RCMIN.EQ.ZERO ) THEN
173: *
174: *        Find the first zero scale factor and return an error code.
175: *
176:          DO 50 I = 1, M
177:             IF( R( I ).EQ.ZERO ) THEN
178:                INFO = I
179:                RETURN
180:             END IF
181:    50    CONTINUE
182:       ELSE
183: *
184: *        Invert the scale factors.
185: *
186:          DO 60 I = 1, M
187:             R( I ) = ONE / MIN( MAX( R( I ), SMLNUM ), BIGNUM )
188:    60    CONTINUE
189: *
190: *        Compute ROWCND = min(R(I)) / max(R(I))
191: *
192:          ROWCND = MAX( RCMIN, SMLNUM ) / MIN( RCMAX, BIGNUM )
193:       END IF
194: *
195: *     Compute column scale factors
196: *
197:       DO 70 J = 1, N
198:          C( J ) = ZERO
199:    70 CONTINUE
200: *
201: *     Find the maximum element in each column,
202: *     assuming the row scaling computed above.
203: *
204:       KD = KU + 1
205:       DO 90 J = 1, N
206:          DO 80 I = MAX( J-KU, 1 ), MIN( J+KL, M )
207:             C( J ) = MAX( C( J ), CABS1( AB( KD+I-J, J ) )*R( I ) )
208:    80    CONTINUE
209:    90 CONTINUE
210: *
211: *     Find the maximum and minimum scale factors.
212: *
213:       RCMIN = BIGNUM
214:       RCMAX = ZERO
215:       DO 100 J = 1, N
216:          RCMIN = MIN( RCMIN, C( J ) )
217:          RCMAX = MAX( RCMAX, C( J ) )
218:   100 CONTINUE
219: *
220:       IF( RCMIN.EQ.ZERO ) THEN
221: *
222: *        Find the first zero scale factor and return an error code.
223: *
224:          DO 110 J = 1, N
225:             IF( C( J ).EQ.ZERO ) THEN
226:                INFO = M + J
227:                RETURN
228:             END IF
229:   110    CONTINUE
230:       ELSE
231: *
232: *        Invert the scale factors.
233: *
234:          DO 120 J = 1, N
235:             C( J ) = ONE / MIN( MAX( C( J ), SMLNUM ), BIGNUM )
236:   120    CONTINUE
237: *
238: *        Compute COLCND = min(C(J)) / max(C(J))
239: *
240:          COLCND = MAX( RCMIN, SMLNUM ) / MIN( RCMAX, BIGNUM )
241:       END IF
242: *
243:       RETURN
244: *
245: *     End of CGBEQU
246: *
247:       END
248: