00001 SUBROUTINE ZGBEQU( M, N, KL, KU, AB, LDAB, R, C, ROWCND, COLCND, 00002 $ AMAX, INFO ) 00003 * 00004 * -- LAPACK routine (version 3.2) -- 00005 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00006 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00007 * November 2006 00008 * 00009 * .. Scalar Arguments .. 00010 INTEGER INFO, KL, KU, LDAB, M, N 00011 DOUBLE PRECISION AMAX, COLCND, ROWCND 00012 * .. 00013 * .. Array Arguments .. 00014 DOUBLE PRECISION C( * ), R( * ) 00015 COMPLEX*16 AB( LDAB, * ) 00016 * .. 00017 * 00018 * Purpose 00019 * ======= 00020 * 00021 * ZGBEQU computes row and column scalings intended to equilibrate an 00022 * M-by-N band matrix A and reduce its condition number. R returns the 00023 * row scale factors and C the column scale factors, chosen to try to 00024 * make the largest element in each row and column of the matrix B with 00025 * elements B(i,j)=R(i)*A(i,j)*C(j) have absolute value 1. 00026 * 00027 * R(i) and C(j) are restricted to be between SMLNUM = smallest safe 00028 * number and BIGNUM = largest safe number. Use of these scaling 00029 * factors is not guaranteed to reduce the condition number of A but 00030 * works well in practice. 00031 * 00032 * Arguments 00033 * ========= 00034 * 00035 * M (input) INTEGER 00036 * The number of rows of the matrix A. M >= 0. 00037 * 00038 * N (input) INTEGER 00039 * The number of columns of the matrix A. N >= 0. 00040 * 00041 * KL (input) INTEGER 00042 * The number of subdiagonals within the band of A. KL >= 0. 00043 * 00044 * KU (input) INTEGER 00045 * The number of superdiagonals within the band of A. KU >= 0. 00046 * 00047 * AB (input) COMPLEX*16 array, dimension (LDAB,N) 00048 * The band matrix A, stored in rows 1 to KL+KU+1. The j-th 00049 * column of A is stored in the j-th column of the array AB as 00050 * follows: 00051 * AB(ku+1+i-j,j) = A(i,j) for max(1,j-ku)<=i<=min(m,j+kl). 00052 * 00053 * LDAB (input) INTEGER 00054 * The leading dimension of the array AB. LDAB >= KL+KU+1. 00055 * 00056 * R (output) DOUBLE PRECISION array, dimension (M) 00057 * If INFO = 0, or INFO > M, R contains the row scale factors 00058 * for A. 00059 * 00060 * C (output) DOUBLE PRECISION array, dimension (N) 00061 * If INFO = 0, C contains the column scale factors for A. 00062 * 00063 * ROWCND (output) DOUBLE PRECISION 00064 * If INFO = 0 or INFO > M, ROWCND contains the ratio of the 00065 * smallest R(i) to the largest R(i). If ROWCND >= 0.1 and 00066 * AMAX is neither too large nor too small, it is not worth 00067 * scaling by R. 00068 * 00069 * COLCND (output) DOUBLE PRECISION 00070 * If INFO = 0, COLCND contains the ratio of the smallest 00071 * C(i) to the largest C(i). If COLCND >= 0.1, it is not 00072 * worth scaling by C. 00073 * 00074 * AMAX (output) DOUBLE PRECISION 00075 * Absolute value of largest matrix element. If AMAX is very 00076 * close to overflow or very close to underflow, the matrix 00077 * should be scaled. 00078 * 00079 * INFO (output) INTEGER 00080 * = 0: successful exit 00081 * < 0: if INFO = -i, the i-th argument had an illegal value 00082 * > 0: if INFO = i, and i is 00083 * <= M: the i-th row of A is exactly zero 00084 * > M: the (i-M)-th column of A is exactly zero 00085 * 00086 * ===================================================================== 00087 * 00088 * .. Parameters .. 00089 DOUBLE PRECISION ONE, ZERO 00090 PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 ) 00091 * .. 00092 * .. Local Scalars .. 00093 INTEGER I, J, KD 00094 DOUBLE PRECISION BIGNUM, RCMAX, RCMIN, SMLNUM 00095 COMPLEX*16 ZDUM 00096 * .. 00097 * .. External Functions .. 00098 DOUBLE PRECISION DLAMCH 00099 EXTERNAL DLAMCH 00100 * .. 00101 * .. External Subroutines .. 00102 EXTERNAL XERBLA 00103 * .. 00104 * .. Intrinsic Functions .. 00105 INTRINSIC ABS, DBLE, DIMAG, MAX, MIN 00106 * .. 00107 * .. Statement Functions .. 00108 DOUBLE PRECISION CABS1 00109 * .. 00110 * .. Statement Function definitions .. 00111 CABS1( ZDUM ) = ABS( DBLE( ZDUM ) ) + ABS( DIMAG( ZDUM ) ) 00112 * .. 00113 * .. Executable Statements .. 00114 * 00115 * Test the input parameters 00116 * 00117 INFO = 0 00118 IF( M.LT.0 ) THEN 00119 INFO = -1 00120 ELSE IF( N.LT.0 ) THEN 00121 INFO = -2 00122 ELSE IF( KL.LT.0 ) THEN 00123 INFO = -3 00124 ELSE IF( KU.LT.0 ) THEN 00125 INFO = -4 00126 ELSE IF( LDAB.LT.KL+KU+1 ) THEN 00127 INFO = -6 00128 END IF 00129 IF( INFO.NE.0 ) THEN 00130 CALL XERBLA( 'ZGBEQU', -INFO ) 00131 RETURN 00132 END IF 00133 * 00134 * Quick return if possible 00135 * 00136 IF( M.EQ.0 .OR. N.EQ.0 ) THEN 00137 ROWCND = ONE 00138 COLCND = ONE 00139 AMAX = ZERO 00140 RETURN 00141 END IF 00142 * 00143 * Get machine constants. 00144 * 00145 SMLNUM = DLAMCH( 'S' ) 00146 BIGNUM = ONE / SMLNUM 00147 * 00148 * Compute row scale factors. 00149 * 00150 DO 10 I = 1, M 00151 R( I ) = ZERO 00152 10 CONTINUE 00153 * 00154 * Find the maximum element in each row. 00155 * 00156 KD = KU + 1 00157 DO 30 J = 1, N 00158 DO 20 I = MAX( J-KU, 1 ), MIN( J+KL, M ) 00159 R( I ) = MAX( R( I ), CABS1( AB( KD+I-J, J ) ) ) 00160 20 CONTINUE 00161 30 CONTINUE 00162 * 00163 * Find the maximum and minimum scale factors. 00164 * 00165 RCMIN = BIGNUM 00166 RCMAX = ZERO 00167 DO 40 I = 1, M 00168 RCMAX = MAX( RCMAX, R( I ) ) 00169 RCMIN = MIN( RCMIN, R( I ) ) 00170 40 CONTINUE 00171 AMAX = RCMAX 00172 * 00173 IF( RCMIN.EQ.ZERO ) THEN 00174 * 00175 * Find the first zero scale factor and return an error code. 00176 * 00177 DO 50 I = 1, M 00178 IF( R( I ).EQ.ZERO ) THEN 00179 INFO = I 00180 RETURN 00181 END IF 00182 50 CONTINUE 00183 ELSE 00184 * 00185 * Invert the scale factors. 00186 * 00187 DO 60 I = 1, M 00188 R( I ) = ONE / MIN( MAX( R( I ), SMLNUM ), BIGNUM ) 00189 60 CONTINUE 00190 * 00191 * Compute ROWCND = min(R(I)) / max(R(I)) 00192 * 00193 ROWCND = MAX( RCMIN, SMLNUM ) / MIN( RCMAX, BIGNUM ) 00194 END IF 00195 * 00196 * Compute column scale factors 00197 * 00198 DO 70 J = 1, N 00199 C( J ) = ZERO 00200 70 CONTINUE 00201 * 00202 * Find the maximum element in each column, 00203 * assuming the row scaling computed above. 00204 * 00205 KD = KU + 1 00206 DO 90 J = 1, N 00207 DO 80 I = MAX( J-KU, 1 ), MIN( J+KL, M ) 00208 C( J ) = MAX( C( J ), CABS1( AB( KD+I-J, J ) )*R( I ) ) 00209 80 CONTINUE 00210 90 CONTINUE 00211 * 00212 * Find the maximum and minimum scale factors. 00213 * 00214 RCMIN = BIGNUM 00215 RCMAX = ZERO 00216 DO 100 J = 1, N 00217 RCMIN = MIN( RCMIN, C( J ) ) 00218 RCMAX = MAX( RCMAX, C( J ) ) 00219 100 CONTINUE 00220 * 00221 IF( RCMIN.EQ.ZERO ) THEN 00222 * 00223 * Find the first zero scale factor and return an error code. 00224 * 00225 DO 110 J = 1, N 00226 IF( C( J ).EQ.ZERO ) THEN 00227 INFO = M + J 00228 RETURN 00229 END IF 00230 110 CONTINUE 00231 ELSE 00232 * 00233 * Invert the scale factors. 00234 * 00235 DO 120 J = 1, N 00236 C( J ) = ONE / MIN( MAX( C( J ), SMLNUM ), BIGNUM ) 00237 120 CONTINUE 00238 * 00239 * Compute COLCND = min(C(J)) / max(C(J)) 00240 * 00241 COLCND = MAX( RCMIN, SMLNUM ) / MIN( RCMAX, BIGNUM ) 00242 END IF 00243 * 00244 RETURN 00245 * 00246 * End of ZGBEQU 00247 * 00248 END
1.7.2