001:       SUBROUTINE ZGEEQU( M, N, A, LDA, R, C, ROWCND, COLCND, AMAX,
002:      $                   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, LDA, M, N
010:       DOUBLE PRECISION   AMAX, COLCND, ROWCND
011: *     ..
012: *     .. Array Arguments ..
013:       DOUBLE PRECISION   C( * ), R( * )
014:       COMPLEX*16         A( LDA, * )
015: *     ..
016: *
017: *  Purpose
018: *  =======
019: *
020: *  ZGEEQU computes row and column scalings intended to equilibrate an
021: *  M-by-N matrix A and reduce its condition number.  R returns the row
022: *  scale factors and C the column scale factors, chosen to try to make
023: *  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: *  A       (input) COMPLEX*16 array, dimension (LDA,N)
041: *          The M-by-N matrix whose equilibration factors are
042: *          to be computed.
043: *
044: *  LDA     (input) INTEGER
045: *          The leading dimension of the array A.  LDA >= max(1,M).
046: *
047: *  R       (output) DOUBLE PRECISION array, dimension (M)
048: *          If INFO = 0 or INFO > M, R contains the row scale factors
049: *          for A.
050: *
051: *  C       (output) DOUBLE PRECISION array, dimension (N)
052: *          If INFO = 0,  C contains the column scale factors for A.
053: *
054: *  ROWCND  (output) DOUBLE PRECISION
055: *          If INFO = 0 or INFO > M, ROWCND contains the ratio of the
056: *          smallest R(i) to the largest R(i).  If ROWCND >= 0.1 and
057: *          AMAX is neither too large nor too small, it is not worth
058: *          scaling by R.
059: *
060: *  COLCND  (output) DOUBLE PRECISION
061: *          If INFO = 0, COLCND contains the ratio of the smallest
062: *          C(i) to the largest C(i).  If COLCND >= 0.1, it is not
063: *          worth scaling by C.
064: *
065: *  AMAX    (output) DOUBLE PRECISION
066: *          Absolute value of largest matrix element.  If AMAX is very
067: *          close to overflow or very close to underflow, the matrix
068: *          should be scaled.
069: *
070: *  INFO    (output) INTEGER
071: *          = 0:  successful exit
072: *          < 0:  if INFO = -i, the i-th argument had an illegal value
073: *          > 0:  if INFO = i,  and i is
074: *                <= M:  the i-th row of A is exactly zero
075: *                >  M:  the (i-M)-th column of A is exactly zero
076: *
077: *  =====================================================================
078: *
079: *     .. Parameters ..
080:       DOUBLE PRECISION   ONE, ZERO
081:       PARAMETER          ( ONE = 1.0D+0, ZERO = 0.0D+0 )
082: *     ..
083: *     .. Local Scalars ..
084:       INTEGER            I, J
085:       DOUBLE PRECISION   BIGNUM, RCMAX, RCMIN, SMLNUM
086:       COMPLEX*16         ZDUM
087: *     ..
088: *     .. External Functions ..
089:       DOUBLE PRECISION   DLAMCH
090:       EXTERNAL           DLAMCH
091: *     ..
092: *     .. External Subroutines ..
093:       EXTERNAL           XERBLA
094: *     ..
095: *     .. Intrinsic Functions ..
096:       INTRINSIC          ABS, DBLE, DIMAG, MAX, MIN
097: *     ..
098: *     .. Statement Functions ..
099:       DOUBLE PRECISION   CABS1
100: *     ..
101: *     .. Statement Function definitions ..
102:       CABS1( ZDUM ) = ABS( DBLE( ZDUM ) ) + ABS( DIMAG( ZDUM ) )
103: *     ..
104: *     .. Executable Statements ..
105: *
106: *     Test the input parameters.
107: *
108:       INFO = 0
109:       IF( M.LT.0 ) THEN
110:          INFO = -1
111:       ELSE IF( N.LT.0 ) THEN
112:          INFO = -2
113:       ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
114:          INFO = -4
115:       END IF
116:       IF( INFO.NE.0 ) THEN
117:          CALL XERBLA( 'ZGEEQU', -INFO )
118:          RETURN
119:       END IF
120: *
121: *     Quick return if possible
122: *
123:       IF( M.EQ.0 .OR. N.EQ.0 ) THEN
124:          ROWCND = ONE
125:          COLCND = ONE
126:          AMAX = ZERO
127:          RETURN
128:       END IF
129: *
130: *     Get machine constants.
131: *
132:       SMLNUM = DLAMCH( 'S' )
133:       BIGNUM = ONE / SMLNUM
134: *
135: *     Compute row scale factors.
136: *
137:       DO 10 I = 1, M
138:          R( I ) = ZERO
139:    10 CONTINUE
140: *
141: *     Find the maximum element in each row.
142: *
143:       DO 30 J = 1, N
144:          DO 20 I = 1, M
145:             R( I ) = MAX( R( I ), CABS1( A( I, J ) ) )
146:    20    CONTINUE
147:    30 CONTINUE
148: *
149: *     Find the maximum and minimum scale factors.
150: *
151:       RCMIN = BIGNUM
152:       RCMAX = ZERO
153:       DO 40 I = 1, M
154:          RCMAX = MAX( RCMAX, R( I ) )
155:          RCMIN = MIN( RCMIN, R( I ) )
156:    40 CONTINUE
157:       AMAX = RCMAX
158: *
159:       IF( RCMIN.EQ.ZERO ) THEN
160: *
161: *        Find the first zero scale factor and return an error code.
162: *
163:          DO 50 I = 1, M
164:             IF( R( I ).EQ.ZERO ) THEN
165:                INFO = I
166:                RETURN
167:             END IF
168:    50    CONTINUE
169:       ELSE
170: *
171: *        Invert the scale factors.
172: *
173:          DO 60 I = 1, M
174:             R( I ) = ONE / MIN( MAX( R( I ), SMLNUM ), BIGNUM )
175:    60    CONTINUE
176: *
177: *        Compute ROWCND = min(R(I)) / max(R(I))
178: *
179:          ROWCND = MAX( RCMIN, SMLNUM ) / MIN( RCMAX, BIGNUM )
180:       END IF
181: *
182: *     Compute column scale factors
183: *
184:       DO 70 J = 1, N
185:          C( J ) = ZERO
186:    70 CONTINUE
187: *
188: *     Find the maximum element in each column,
189: *     assuming the row scaling computed above.
190: *
191:       DO 90 J = 1, N
192:          DO 80 I = 1, M
193:             C( J ) = MAX( C( J ), CABS1( A( I, J ) )*R( I ) )
194:    80    CONTINUE
195:    90 CONTINUE
196: *
197: *     Find the maximum and minimum scale factors.
198: *
199:       RCMIN = BIGNUM
200:       RCMAX = ZERO
201:       DO 100 J = 1, N
202:          RCMIN = MIN( RCMIN, C( J ) )
203:          RCMAX = MAX( RCMAX, C( J ) )
204:   100 CONTINUE
205: *
206:       IF( RCMIN.EQ.ZERO ) THEN
207: *
208: *        Find the first zero scale factor and return an error code.
209: *
210:          DO 110 J = 1, N
211:             IF( C( J ).EQ.ZERO ) THEN
212:                INFO = M + J
213:                RETURN
214:             END IF
215:   110    CONTINUE
216:       ELSE
217: *
218: *        Invert the scale factors.
219: *
220:          DO 120 J = 1, N
221:             C( J ) = ONE / MIN( MAX( C( J ), SMLNUM ), BIGNUM )
222:   120    CONTINUE
223: *
224: *        Compute COLCND = min(C(J)) / max(C(J))
225: *
226:          COLCND = MAX( RCMIN, SMLNUM ) / MIN( RCMAX, BIGNUM )
227:       END IF
228: *
229:       RETURN
230: *
231: *     End of ZGEEQU
232: *
233:       END
234: