1 DOUBLE PRECISION FUNCTION ZLA_GERCOND_X( TRANS, N, A, LDA, AF,
2 $ LDAF, IPIV, X, INFO,
3 $ WORK, RWORK )
4 *
5 * -- LAPACK routine (version 3.2.1) --
6 * -- Contributed by James Demmel, Deaglan Halligan, Yozo Hida and --
7 * -- Jason Riedy of Univ. of California Berkeley. --
8 * -- April 2009 --
9 *
10 * -- LAPACK is a software package provided by Univ. of Tennessee, --
11 * -- Univ. of California Berkeley and NAG Ltd. --
12 *
13 IMPLICIT NONE
14 * ..
15 * .. Scalar Arguments ..
16 CHARACTER TRANS
17 INTEGER N, LDA, LDAF, INFO
18 * ..
19 * .. Array Arguments ..
20 INTEGER IPIV( * )
21 COMPLEX*16 A( LDA, * ), AF( LDAF, * ), WORK( * ), X( * )
22 DOUBLE PRECISION RWORK( * )
23 * ..
24 *
25 * Purpose
26 * =======
27 *
28 * ZLA_GERCOND_X computes the infinity norm condition number of
29 * op(A) * diag(X) where X is a COMPLEX*16 vector.
30 *
31 * Arguments
32 * =========
33 *
34 * TRANS (input) CHARACTER*1
35 * Specifies the form of the system of equations:
36 * = 'N': A * X = B (No transpose)
37 * = 'T': A**T * X = B (Transpose)
38 * = 'C': A**H * X = B (Conjugate Transpose = Transpose)
39 *
40 * N (input) INTEGER
41 * The number of linear equations, i.e., the order of the
42 * matrix A. N >= 0.
43 *
44 * A (input) COMPLEX*16 array, dimension (LDA,N)
45 * On entry, the N-by-N matrix A.
46 *
47 * LDA (input) INTEGER
48 * The leading dimension of the array A. LDA >= max(1,N).
49 *
50 * AF (input) COMPLEX*16 array, dimension (LDAF,N)
51 * The factors L and U from the factorization
52 * A = P*L*U as computed by ZGETRF.
53 *
54 * LDAF (input) INTEGER
55 * The leading dimension of the array AF. LDAF >= max(1,N).
56 *
57 * IPIV (input) INTEGER array, dimension (N)
58 * The pivot indices from the factorization A = P*L*U
59 * as computed by ZGETRF; row i of the matrix was interchanged
60 * with row IPIV(i).
61 *
62 * X (input) COMPLEX*16 array, dimension (N)
63 * The vector X in the formula op(A) * diag(X).
64 *
65 * INFO (output) INTEGER
66 * = 0: Successful exit.
67 * i > 0: The ith argument is invalid.
68 *
69 * WORK (input) COMPLEX*16 array, dimension (2*N).
70 * Workspace.
71 *
72 * RWORK (input) DOUBLE PRECISION array, dimension (N).
73 * Workspace.
74 *
75 * =====================================================================
76 *
77 * .. Local Scalars ..
78 LOGICAL NOTRANS
79 INTEGER KASE
80 DOUBLE PRECISION AINVNM, ANORM, TMP
81 INTEGER I, J
82 COMPLEX*16 ZDUM
83 * ..
84 * .. Local Arrays ..
85 INTEGER ISAVE( 3 )
86 * ..
87 * .. External Functions ..
88 LOGICAL LSAME
89 EXTERNAL LSAME
90 * ..
91 * .. External Subroutines ..
92 EXTERNAL ZLACN2, ZGETRS, XERBLA
93 * ..
94 * .. Intrinsic Functions ..
95 INTRINSIC ABS, MAX, REAL, DIMAG
96 * ..
97 * .. Statement Functions ..
98 DOUBLE PRECISION CABS1
99 * ..
100 * .. Statement Function Definitions ..
101 CABS1( ZDUM ) = ABS( DBLE( ZDUM ) ) + ABS( DIMAG( ZDUM ) )
102 * ..
103 * .. Executable Statements ..
104 *
105 ZLA_GERCOND_X = 0.0D+0
106 *
107 INFO = 0
108 NOTRANS = LSAME( TRANS, 'N' )
109 IF ( .NOT. NOTRANS .AND. .NOT. LSAME( TRANS, 'T' ) .AND. .NOT.
110 $ LSAME( TRANS, 'C' ) ) THEN
111 INFO = -1
112 ELSE IF( N.LT.0 ) THEN
113 INFO = -2
114 END IF
115 IF( INFO.NE.0 ) THEN
116 CALL XERBLA( 'ZLA_GERCOND_X', -INFO )
117 RETURN
118 END IF
119 *
120 * Compute norm of op(A)*op2(C).
121 *
122 ANORM = 0.0D+0
123 IF ( NOTRANS ) THEN
124 DO I = 1, N
125 TMP = 0.0D+0
126 DO J = 1, N
127 TMP = TMP + CABS1( A( I, J ) * X( J ) )
128 END DO
129 RWORK( I ) = TMP
130 ANORM = MAX( ANORM, TMP )
131 END DO
132 ELSE
133 DO I = 1, N
134 TMP = 0.0D+0
135 DO J = 1, N
136 TMP = TMP + CABS1( A( J, I ) * X( J ) )
137 END DO
138 RWORK( I ) = TMP
139 ANORM = MAX( ANORM, TMP )
140 END DO
141 END IF
142 *
143 * Quick return if possible.
144 *
145 IF( N.EQ.0 ) THEN
146 ZLA_GERCOND_X = 1.0D+0
147 RETURN
148 ELSE IF( ANORM .EQ. 0.0D+0 ) THEN
149 RETURN
150 END IF
151 *
152 * Estimate the norm of inv(op(A)).
153 *
154 AINVNM = 0.0D+0
155 *
156 KASE = 0
157 10 CONTINUE
158 CALL ZLACN2( N, WORK( N+1 ), WORK, AINVNM, KASE, ISAVE )
159 IF( KASE.NE.0 ) THEN
160 IF( KASE.EQ.2 ) THEN
161 * Multiply by R.
162 DO I = 1, N
163 WORK( I ) = WORK( I ) * RWORK( I )
164 END DO
165 *
166 IF ( NOTRANS ) THEN
167 CALL ZGETRS( 'No transpose', N, 1, AF, LDAF, IPIV,
168 $ WORK, N, INFO )
169 ELSE
170 CALL ZGETRS( 'Conjugate transpose', N, 1, AF, LDAF, IPIV,
171 $ WORK, N, INFO )
172 ENDIF
173 *
174 * Multiply by inv(X).
175 *
176 DO I = 1, N
177 WORK( I ) = WORK( I ) / X( I )
178 END DO
179 ELSE
180 *
181 * Multiply by inv(X**H).
182 *
183 DO I = 1, N
184 WORK( I ) = WORK( I ) / X( I )
185 END DO
186 *
187 IF ( NOTRANS ) THEN
188 CALL ZGETRS( 'Conjugate transpose', N, 1, AF, LDAF, IPIV,
189 $ WORK, N, INFO )
190 ELSE
191 CALL ZGETRS( 'No transpose', N, 1, AF, LDAF, IPIV,
192 $ WORK, N, INFO )
193 END IF
194 *
195 * Multiply by R.
196 *
197 DO I = 1, N
198 WORK( I ) = WORK( I ) * RWORK( I )
199 END DO
200 END IF
201 GO TO 10
202 END IF
203 *
204 * Compute the estimate of the reciprocal condition number.
205 *
206 IF( AINVNM .NE. 0.0D+0 )
207 $ ZLA_GERCOND_X = 1.0D+0 / AINVNM
208 *
209 RETURN
210 *
211 END
2 $ LDAF, IPIV, X, INFO,
3 $ WORK, RWORK )
4 *
5 * -- LAPACK routine (version 3.2.1) --
6 * -- Contributed by James Demmel, Deaglan Halligan, Yozo Hida and --
7 * -- Jason Riedy of Univ. of California Berkeley. --
8 * -- April 2009 --
9 *
10 * -- LAPACK is a software package provided by Univ. of Tennessee, --
11 * -- Univ. of California Berkeley and NAG Ltd. --
12 *
13 IMPLICIT NONE
14 * ..
15 * .. Scalar Arguments ..
16 CHARACTER TRANS
17 INTEGER N, LDA, LDAF, INFO
18 * ..
19 * .. Array Arguments ..
20 INTEGER IPIV( * )
21 COMPLEX*16 A( LDA, * ), AF( LDAF, * ), WORK( * ), X( * )
22 DOUBLE PRECISION RWORK( * )
23 * ..
24 *
25 * Purpose
26 * =======
27 *
28 * ZLA_GERCOND_X computes the infinity norm condition number of
29 * op(A) * diag(X) where X is a COMPLEX*16 vector.
30 *
31 * Arguments
32 * =========
33 *
34 * TRANS (input) CHARACTER*1
35 * Specifies the form of the system of equations:
36 * = 'N': A * X = B (No transpose)
37 * = 'T': A**T * X = B (Transpose)
38 * = 'C': A**H * X = B (Conjugate Transpose = Transpose)
39 *
40 * N (input) INTEGER
41 * The number of linear equations, i.e., the order of the
42 * matrix A. N >= 0.
43 *
44 * A (input) COMPLEX*16 array, dimension (LDA,N)
45 * On entry, the N-by-N matrix A.
46 *
47 * LDA (input) INTEGER
48 * The leading dimension of the array A. LDA >= max(1,N).
49 *
50 * AF (input) COMPLEX*16 array, dimension (LDAF,N)
51 * The factors L and U from the factorization
52 * A = P*L*U as computed by ZGETRF.
53 *
54 * LDAF (input) INTEGER
55 * The leading dimension of the array AF. LDAF >= max(1,N).
56 *
57 * IPIV (input) INTEGER array, dimension (N)
58 * The pivot indices from the factorization A = P*L*U
59 * as computed by ZGETRF; row i of the matrix was interchanged
60 * with row IPIV(i).
61 *
62 * X (input) COMPLEX*16 array, dimension (N)
63 * The vector X in the formula op(A) * diag(X).
64 *
65 * INFO (output) INTEGER
66 * = 0: Successful exit.
67 * i > 0: The ith argument is invalid.
68 *
69 * WORK (input) COMPLEX*16 array, dimension (2*N).
70 * Workspace.
71 *
72 * RWORK (input) DOUBLE PRECISION array, dimension (N).
73 * Workspace.
74 *
75 * =====================================================================
76 *
77 * .. Local Scalars ..
78 LOGICAL NOTRANS
79 INTEGER KASE
80 DOUBLE PRECISION AINVNM, ANORM, TMP
81 INTEGER I, J
82 COMPLEX*16 ZDUM
83 * ..
84 * .. Local Arrays ..
85 INTEGER ISAVE( 3 )
86 * ..
87 * .. External Functions ..
88 LOGICAL LSAME
89 EXTERNAL LSAME
90 * ..
91 * .. External Subroutines ..
92 EXTERNAL ZLACN2, ZGETRS, XERBLA
93 * ..
94 * .. Intrinsic Functions ..
95 INTRINSIC ABS, MAX, REAL, DIMAG
96 * ..
97 * .. Statement Functions ..
98 DOUBLE PRECISION CABS1
99 * ..
100 * .. Statement Function Definitions ..
101 CABS1( ZDUM ) = ABS( DBLE( ZDUM ) ) + ABS( DIMAG( ZDUM ) )
102 * ..
103 * .. Executable Statements ..
104 *
105 ZLA_GERCOND_X = 0.0D+0
106 *
107 INFO = 0
108 NOTRANS = LSAME( TRANS, 'N' )
109 IF ( .NOT. NOTRANS .AND. .NOT. LSAME( TRANS, 'T' ) .AND. .NOT.
110 $ LSAME( TRANS, 'C' ) ) THEN
111 INFO = -1
112 ELSE IF( N.LT.0 ) THEN
113 INFO = -2
114 END IF
115 IF( INFO.NE.0 ) THEN
116 CALL XERBLA( 'ZLA_GERCOND_X', -INFO )
117 RETURN
118 END IF
119 *
120 * Compute norm of op(A)*op2(C).
121 *
122 ANORM = 0.0D+0
123 IF ( NOTRANS ) THEN
124 DO I = 1, N
125 TMP = 0.0D+0
126 DO J = 1, N
127 TMP = TMP + CABS1( A( I, J ) * X( J ) )
128 END DO
129 RWORK( I ) = TMP
130 ANORM = MAX( ANORM, TMP )
131 END DO
132 ELSE
133 DO I = 1, N
134 TMP = 0.0D+0
135 DO J = 1, N
136 TMP = TMP + CABS1( A( J, I ) * X( J ) )
137 END DO
138 RWORK( I ) = TMP
139 ANORM = MAX( ANORM, TMP )
140 END DO
141 END IF
142 *
143 * Quick return if possible.
144 *
145 IF( N.EQ.0 ) THEN
146 ZLA_GERCOND_X = 1.0D+0
147 RETURN
148 ELSE IF( ANORM .EQ. 0.0D+0 ) THEN
149 RETURN
150 END IF
151 *
152 * Estimate the norm of inv(op(A)).
153 *
154 AINVNM = 0.0D+0
155 *
156 KASE = 0
157 10 CONTINUE
158 CALL ZLACN2( N, WORK( N+1 ), WORK, AINVNM, KASE, ISAVE )
159 IF( KASE.NE.0 ) THEN
160 IF( KASE.EQ.2 ) THEN
161 * Multiply by R.
162 DO I = 1, N
163 WORK( I ) = WORK( I ) * RWORK( I )
164 END DO
165 *
166 IF ( NOTRANS ) THEN
167 CALL ZGETRS( 'No transpose', N, 1, AF, LDAF, IPIV,
168 $ WORK, N, INFO )
169 ELSE
170 CALL ZGETRS( 'Conjugate transpose', N, 1, AF, LDAF, IPIV,
171 $ WORK, N, INFO )
172 ENDIF
173 *
174 * Multiply by inv(X).
175 *
176 DO I = 1, N
177 WORK( I ) = WORK( I ) / X( I )
178 END DO
179 ELSE
180 *
181 * Multiply by inv(X**H).
182 *
183 DO I = 1, N
184 WORK( I ) = WORK( I ) / X( I )
185 END DO
186 *
187 IF ( NOTRANS ) THEN
188 CALL ZGETRS( 'Conjugate transpose', N, 1, AF, LDAF, IPIV,
189 $ WORK, N, INFO )
190 ELSE
191 CALL ZGETRS( 'No transpose', N, 1, AF, LDAF, IPIV,
192 $ WORK, N, INFO )
193 END IF
194 *
195 * Multiply by R.
196 *
197 DO I = 1, N
198 WORK( I ) = WORK( I ) * RWORK( I )
199 END DO
200 END IF
201 GO TO 10
202 END IF
203 *
204 * Compute the estimate of the reciprocal condition number.
205 *
206 IF( AINVNM .NE. 0.0D+0 )
207 $ ZLA_GERCOND_X = 1.0D+0 / AINVNM
208 *
209 RETURN
210 *
211 END