Annotation of rpl/lapack/lapack/zgeequ.f, revision 1.2
1.1 bertrand 1: SUBROUTINE ZGEEQU( M, N, A, LDA, R, C, ROWCND, COLCND, AMAX,
2: $ INFO )
3: *
4: * -- LAPACK routine (version 3.2) --
5: * -- LAPACK is a software package provided by Univ. of Tennessee, --
6: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
7: * November 2006
8: *
9: * .. Scalar Arguments ..
10: INTEGER INFO, LDA, M, N
11: DOUBLE PRECISION AMAX, COLCND, ROWCND
12: * ..
13: * .. Array Arguments ..
14: DOUBLE PRECISION C( * ), R( * )
15: COMPLEX*16 A( LDA, * )
16: * ..
17: *
18: * Purpose
19: * =======
20: *
21: * ZGEEQU computes row and column scalings intended to equilibrate an
22: * M-by-N matrix A and reduce its condition number. R returns the row
23: * scale factors and C the column scale factors, chosen to try to make
24: * the largest element in each row and column of the matrix B with
25: * elements B(i,j)=R(i)*A(i,j)*C(j) have absolute value 1.
26: *
27: * R(i) and C(j) are restricted to be between SMLNUM = smallest safe
28: * number and BIGNUM = largest safe number. Use of these scaling
29: * factors is not guaranteed to reduce the condition number of A but
30: * works well in practice.
31: *
32: * Arguments
33: * =========
34: *
35: * M (input) INTEGER
36: * The number of rows of the matrix A. M >= 0.
37: *
38: * N (input) INTEGER
39: * The number of columns of the matrix A. N >= 0.
40: *
41: * A (input) COMPLEX*16 array, dimension (LDA,N)
42: * The M-by-N matrix whose equilibration factors are
43: * to be computed.
44: *
45: * LDA (input) INTEGER
46: * The leading dimension of the array A. LDA >= max(1,M).
47: *
48: * R (output) DOUBLE PRECISION array, dimension (M)
49: * If INFO = 0 or INFO > M, R contains the row scale factors
50: * for A.
51: *
52: * C (output) DOUBLE PRECISION array, dimension (N)
53: * If INFO = 0, C contains the column scale factors for A.
54: *
55: * ROWCND (output) DOUBLE PRECISION
56: * If INFO = 0 or INFO > M, ROWCND contains the ratio of the
57: * smallest R(i) to the largest R(i). If ROWCND >= 0.1 and
58: * AMAX is neither too large nor too small, it is not worth
59: * scaling by R.
60: *
61: * COLCND (output) DOUBLE PRECISION
62: * If INFO = 0, COLCND contains the ratio of the smallest
63: * C(i) to the largest C(i). If COLCND >= 0.1, it is not
64: * worth scaling by C.
65: *
66: * AMAX (output) DOUBLE PRECISION
67: * Absolute value of largest matrix element. If AMAX is very
68: * close to overflow or very close to underflow, the matrix
69: * should be scaled.
70: *
71: * INFO (output) INTEGER
72: * = 0: successful exit
73: * < 0: if INFO = -i, the i-th argument had an illegal value
74: * > 0: if INFO = i, and i is
75: * <= M: the i-th row of A is exactly zero
76: * > M: the (i-M)-th column of A is exactly zero
77: *
78: * =====================================================================
79: *
80: * .. Parameters ..
81: DOUBLE PRECISION ONE, ZERO
82: PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 )
83: * ..
84: * .. Local Scalars ..
85: INTEGER I, J
86: DOUBLE PRECISION BIGNUM, RCMAX, RCMIN, SMLNUM
87: COMPLEX*16 ZDUM
88: * ..
89: * .. External Functions ..
90: DOUBLE PRECISION DLAMCH
91: EXTERNAL DLAMCH
92: * ..
93: * .. External Subroutines ..
94: EXTERNAL XERBLA
95: * ..
96: * .. Intrinsic Functions ..
97: INTRINSIC ABS, DBLE, DIMAG, MAX, MIN
98: * ..
99: * .. Statement Functions ..
100: DOUBLE PRECISION CABS1
101: * ..
102: * .. Statement Function definitions ..
103: CABS1( ZDUM ) = ABS( DBLE( ZDUM ) ) + ABS( DIMAG( ZDUM ) )
104: * ..
105: * .. Executable Statements ..
106: *
107: * Test the input parameters.
108: *
109: INFO = 0
110: IF( M.LT.0 ) THEN
111: INFO = -1
112: ELSE IF( N.LT.0 ) THEN
113: INFO = -2
114: ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
115: INFO = -4
116: END IF
117: IF( INFO.NE.0 ) THEN
118: CALL XERBLA( 'ZGEEQU', -INFO )
119: RETURN
120: END IF
121: *
122: * Quick return if possible
123: *
124: IF( M.EQ.0 .OR. N.EQ.0 ) THEN
125: ROWCND = ONE
126: COLCND = ONE
127: AMAX = ZERO
128: RETURN
129: END IF
130: *
131: * Get machine constants.
132: *
133: SMLNUM = DLAMCH( 'S' )
134: BIGNUM = ONE / SMLNUM
135: *
136: * Compute row scale factors.
137: *
138: DO 10 I = 1, M
139: R( I ) = ZERO
140: 10 CONTINUE
141: *
142: * Find the maximum element in each row.
143: *
144: DO 30 J = 1, N
145: DO 20 I = 1, M
146: R( I ) = MAX( R( I ), CABS1( A( I, J ) ) )
147: 20 CONTINUE
148: 30 CONTINUE
149: *
150: * Find the maximum and minimum scale factors.
151: *
152: RCMIN = BIGNUM
153: RCMAX = ZERO
154: DO 40 I = 1, M
155: RCMAX = MAX( RCMAX, R( I ) )
156: RCMIN = MIN( RCMIN, R( I ) )
157: 40 CONTINUE
158: AMAX = RCMAX
159: *
160: IF( RCMIN.EQ.ZERO ) THEN
161: *
162: * Find the first zero scale factor and return an error code.
163: *
164: DO 50 I = 1, M
165: IF( R( I ).EQ.ZERO ) THEN
166: INFO = I
167: RETURN
168: END IF
169: 50 CONTINUE
170: ELSE
171: *
172: * Invert the scale factors.
173: *
174: DO 60 I = 1, M
175: R( I ) = ONE / MIN( MAX( R( I ), SMLNUM ), BIGNUM )
176: 60 CONTINUE
177: *
178: * Compute ROWCND = min(R(I)) / max(R(I))
179: *
180: ROWCND = MAX( RCMIN, SMLNUM ) / MIN( RCMAX, BIGNUM )
181: END IF
182: *
183: * Compute column scale factors
184: *
185: DO 70 J = 1, N
186: C( J ) = ZERO
187: 70 CONTINUE
188: *
189: * Find the maximum element in each column,
190: * assuming the row scaling computed above.
191: *
192: DO 90 J = 1, N
193: DO 80 I = 1, M
194: C( J ) = MAX( C( J ), CABS1( A( I, J ) )*R( I ) )
195: 80 CONTINUE
196: 90 CONTINUE
197: *
198: * Find the maximum and minimum scale factors.
199: *
200: RCMIN = BIGNUM
201: RCMAX = ZERO
202: DO 100 J = 1, N
203: RCMIN = MIN( RCMIN, C( J ) )
204: RCMAX = MAX( RCMAX, C( J ) )
205: 100 CONTINUE
206: *
207: IF( RCMIN.EQ.ZERO ) THEN
208: *
209: * Find the first zero scale factor and return an error code.
210: *
211: DO 110 J = 1, N
212: IF( C( J ).EQ.ZERO ) THEN
213: INFO = M + J
214: RETURN
215: END IF
216: 110 CONTINUE
217: ELSE
218: *
219: * Invert the scale factors.
220: *
221: DO 120 J = 1, N
222: C( J ) = ONE / MIN( MAX( C( J ), SMLNUM ), BIGNUM )
223: 120 CONTINUE
224: *
225: * Compute COLCND = min(C(J)) / max(C(J))
226: *
227: COLCND = MAX( RCMIN, SMLNUM ) / MIN( RCMAX, BIGNUM )
228: END IF
229: *
230: RETURN
231: *
232: * End of ZGEEQU
233: *
234: END
CVSweb interface <joel.bertrand@systella.fr>