001:       SUBROUTINE SGEEQU( M, N, A, LDA, R, C, ROWCND, COLCND, AMAX,
002:      $                   INFO )
003: *
004: *  -- LAPACK routine (version 3.2) --
005: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
006: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
007: *     November 2006
008: *
009: *     .. Scalar Arguments ..
010:       INTEGER            INFO, LDA, M, N
011:       REAL               AMAX, COLCND, ROWCND
012: *     ..
013: *     .. Array Arguments ..
014:       REAL               A( LDA, * ), C( * ), R( * )
015: *     ..
016: *
017: *  Purpose
018: *  =======
019: *
020: *  SGEEQU 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) REAL 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) REAL array, dimension (M)
048: *          If INFO = 0 or INFO > M, R contains the row scale factors
049: *          for A.
050: *
051: *  C       (output) REAL array, dimension (N)
052: *          If INFO = 0,  C contains the column scale factors for A.
053: *
054: *  ROWCND  (output) REAL
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) REAL
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) REAL
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:       REAL               ONE, ZERO
081:       PARAMETER          ( ONE = 1.0E+0, ZERO = 0.0E+0 )
082: *     ..
083: *     .. Local Scalars ..
084:       INTEGER            I, J
085:       REAL               BIGNUM, RCMAX, RCMIN, SMLNUM
086: *     ..
087: *     .. External Functions ..
088:       REAL               SLAMCH
089:       EXTERNAL           SLAMCH
090: *     ..
091: *     .. External Subroutines ..
092:       EXTERNAL           XERBLA
093: *     ..
094: *     .. Intrinsic Functions ..
095:       INTRINSIC          ABS, MAX, MIN
096: *     ..
097: *     .. Executable Statements ..
098: *
099: *     Test the input parameters.
100: *
101:       INFO = 0
102:       IF( M.LT.0 ) THEN
103:          INFO = -1
104:       ELSE IF( N.LT.0 ) THEN
105:          INFO = -2
106:       ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
107:          INFO = -4
108:       END IF
109:       IF( INFO.NE.0 ) THEN
110:          CALL XERBLA( 'SGEEQU', -INFO )
111:          RETURN
112:       END IF
113: *
114: *     Quick return if possible
115: *
116:       IF( M.EQ.0 .OR. N.EQ.0 ) THEN
117:          ROWCND = ONE
118:          COLCND = ONE
119:          AMAX = ZERO
120:          RETURN
121:       END IF
122: *
123: *     Get machine constants.
124: *
125:       SMLNUM = SLAMCH( 'S' )
126:       BIGNUM = ONE / SMLNUM
127: *
128: *     Compute row scale factors.
129: *
130:       DO 10 I = 1, M
131:          R( I ) = ZERO
132:    10 CONTINUE
133: *
134: *     Find the maximum element in each row.
135: *
136:       DO 30 J = 1, N
137:          DO 20 I = 1, M
138:             R( I ) = MAX( R( I ), ABS( A( I, J ) ) )
139:    20    CONTINUE
140:    30 CONTINUE
141: *
142: *     Find the maximum and minimum scale factors.
143: *
144:       RCMIN = BIGNUM
145:       RCMAX = ZERO
146:       DO 40 I = 1, M
147:          RCMAX = MAX( RCMAX, R( I ) )
148:          RCMIN = MIN( RCMIN, R( I ) )
149:    40 CONTINUE
150:       AMAX = RCMAX
151: *
152:       IF( RCMIN.EQ.ZERO ) THEN
153: *
154: *        Find the first zero scale factor and return an error code.
155: *
156:          DO 50 I = 1, M
157:             IF( R( I ).EQ.ZERO ) THEN
158:                INFO = I
159:                RETURN
160:             END IF
161:    50    CONTINUE
162:       ELSE
163: *
164: *        Invert the scale factors.
165: *
166:          DO 60 I = 1, M
167:             R( I ) = ONE / MIN( MAX( R( I ), SMLNUM ), BIGNUM )
168:    60    CONTINUE
169: *
170: *        Compute ROWCND = min(R(I)) / max(R(I))
171: *
172:          ROWCND = MAX( RCMIN, SMLNUM ) / MIN( RCMAX, BIGNUM )
173:       END IF
174: *
175: *     Compute column scale factors
176: *
177:       DO 70 J = 1, N
178:          C( J ) = ZERO
179:    70 CONTINUE
180: *
181: *     Find the maximum element in each column,
182: *     assuming the row scaling computed above.
183: *
184:       DO 90 J = 1, N
185:          DO 80 I = 1, M
186:             C( J ) = MAX( C( J ), ABS( A( I, J ) )*R( I ) )
187:    80    CONTINUE
188:    90 CONTINUE
189: *
190: *     Find the maximum and minimum scale factors.
191: *
192:       RCMIN = BIGNUM
193:       RCMAX = ZERO
194:       DO 100 J = 1, N
195:          RCMIN = MIN( RCMIN, C( J ) )
196:          RCMAX = MAX( RCMAX, C( J ) )
197:   100 CONTINUE
198: *
199:       IF( RCMIN.EQ.ZERO ) THEN
200: *
201: *        Find the first zero scale factor and return an error code.
202: *
203:          DO 110 J = 1, N
204:             IF( C( J ).EQ.ZERO ) THEN
205:                INFO = M + J
206:                RETURN
207:             END IF
208:   110    CONTINUE
209:       ELSE
210: *
211: *        Invert the scale factors.
212: *
213:          DO 120 J = 1, N
214:             C( J ) = ONE / MIN( MAX( C( J ), SMLNUM ), BIGNUM )
215:   120    CONTINUE
216: *
217: *        Compute COLCND = min(C(J)) / max(C(J))
218: *
219:          COLCND = MAX( RCMIN, SMLNUM ) / MIN( RCMAX, BIGNUM )
220:       END IF
221: *
222:       RETURN
223: *
224: *     End of SGEEQU
225: *
226:       END
227: