1: *> \brief \b ZLA_PORCOND_X computes the infinity norm condition number of op(A)*diag(x) for Hermitian positive-definite matrices.
2: *
3: * =========== DOCUMENTATION ===========
4: *
5: * Online html documentation available at
6: * http://www.netlib.org/lapack/explore-html/
7: *
8: *> \htmlonly
9: *> Download ZLA_PORCOND_X + dependencies
10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zla_porcond_x.f">
11: *> [TGZ]</a>
12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zla_porcond_x.f">
13: *> [ZIP]</a>
14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zla_porcond_x.f">
15: *> [TXT]</a>
16: *> \endhtmlonly
17: *
18: * Definition:
19: * ===========
20: *
21: * DOUBLE PRECISION FUNCTION ZLA_PORCOND_X( UPLO, N, A, LDA, AF,
22: * LDAF, X, INFO, WORK,
23: * RWORK )
24: *
25: * .. Scalar Arguments ..
26: * CHARACTER UPLO
27: * INTEGER N, LDA, LDAF, INFO
28: * ..
29: * .. Array Arguments ..
30: * COMPLEX*16 A( LDA, * ), AF( LDAF, * ), WORK( * ), X( * )
31: * DOUBLE PRECISION RWORK( * )
32: * ..
33: *
34: *
35: *> \par Purpose:
36: * =============
37: *>
38: *> \verbatim
39: *>
40: *> ZLA_PORCOND_X Computes the infinity norm condition number of
41: *> op(A) * diag(X) where X is a COMPLEX*16 vector.
42: *> \endverbatim
43: *
44: * Arguments:
45: * ==========
46: *
47: *> \param[in] UPLO
48: *> \verbatim
49: *> UPLO is CHARACTER*1
50: *> = 'U': Upper triangle of A is stored;
51: *> = 'L': Lower triangle of A is stored.
52: *> \endverbatim
53: *>
54: *> \param[in] N
55: *> \verbatim
56: *> N is INTEGER
57: *> The number of linear equations, i.e., the order of the
58: *> matrix A. N >= 0.
59: *> \endverbatim
60: *>
61: *> \param[in] A
62: *> \verbatim
63: *> A is COMPLEX*16 array, dimension (LDA,N)
64: *> On entry, the N-by-N matrix A.
65: *> \endverbatim
66: *>
67: *> \param[in] LDA
68: *> \verbatim
69: *> LDA is INTEGER
70: *> The leading dimension of the array A. LDA >= max(1,N).
71: *> \endverbatim
72: *>
73: *> \param[in] AF
74: *> \verbatim
75: *> AF is COMPLEX*16 array, dimension (LDAF,N)
76: *> The triangular factor U or L from the Cholesky factorization
77: *> A = U**H*U or A = L*L**H, as computed by ZPOTRF.
78: *> \endverbatim
79: *>
80: *> \param[in] LDAF
81: *> \verbatim
82: *> LDAF is INTEGER
83: *> The leading dimension of the array AF. LDAF >= max(1,N).
84: *> \endverbatim
85: *>
86: *> \param[in] X
87: *> \verbatim
88: *> X is COMPLEX*16 array, dimension (N)
89: *> The vector X in the formula op(A) * diag(X).
90: *> \endverbatim
91: *>
92: *> \param[out] INFO
93: *> \verbatim
94: *> INFO is INTEGER
95: *> = 0: Successful exit.
96: *> i > 0: The ith argument is invalid.
97: *> \endverbatim
98: *>
99: *> \param[out] WORK
100: *> \verbatim
101: *> WORK is COMPLEX*16 array, dimension (2*N).
102: *> Workspace.
103: *> \endverbatim
104: *>
105: *> \param[out] RWORK
106: *> \verbatim
107: *> RWORK is DOUBLE PRECISION array, dimension (N).
108: *> Workspace.
109: *> \endverbatim
110: *
111: * Authors:
112: * ========
113: *
114: *> \author Univ. of Tennessee
115: *> \author Univ. of California Berkeley
116: *> \author Univ. of Colorado Denver
117: *> \author NAG Ltd.
118: *
119: *> \ingroup complex16POcomputational
120: *
121: * =====================================================================
122: DOUBLE PRECISION FUNCTION ZLA_PORCOND_X( UPLO, N, A, LDA, AF,
123: $ LDAF, X, INFO, WORK,
124: $ RWORK )
125: *
126: * -- LAPACK computational routine --
127: * -- LAPACK is a software package provided by Univ. of Tennessee, --
128: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
129: *
130: * .. Scalar Arguments ..
131: CHARACTER UPLO
132: INTEGER N, LDA, LDAF, INFO
133: * ..
134: * .. Array Arguments ..
135: COMPLEX*16 A( LDA, * ), AF( LDAF, * ), WORK( * ), X( * )
136: DOUBLE PRECISION RWORK( * )
137: * ..
138: *
139: * =====================================================================
140: *
141: * .. Local Scalars ..
142: INTEGER KASE, I, J
143: DOUBLE PRECISION AINVNM, ANORM, TMP
144: LOGICAL UP, UPPER
145: COMPLEX*16 ZDUM
146: * ..
147: * .. Local Arrays ..
148: INTEGER ISAVE( 3 )
149: * ..
150: * .. External Functions ..
151: LOGICAL LSAME
152: EXTERNAL LSAME
153: * ..
154: * .. External Subroutines ..
155: EXTERNAL ZLACN2, ZPOTRS, XERBLA
156: * ..
157: * .. Intrinsic Functions ..
158: INTRINSIC ABS, MAX, REAL, DIMAG
159: * ..
160: * .. Statement Functions ..
161: DOUBLE PRECISION CABS1
162: * ..
163: * .. Statement Function Definitions ..
164: CABS1( ZDUM ) = ABS( DBLE( ZDUM ) ) + ABS( DIMAG( ZDUM ) )
165: * ..
166: * .. Executable Statements ..
167: *
168: ZLA_PORCOND_X = 0.0D+0
169: *
170: INFO = 0
171: UPPER = LSAME( UPLO, 'U' )
172: IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
173: INFO = -1
174: ELSE IF ( N.LT.0 ) THEN
175: INFO = -2
176: ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
177: INFO = -4
178: ELSE IF( LDAF.LT.MAX( 1, N ) ) THEN
179: INFO = -6
180: END IF
181: IF( INFO.NE.0 ) THEN
182: CALL XERBLA( 'ZLA_PORCOND_X', -INFO )
183: RETURN
184: END IF
185: UP = .FALSE.
186: IF ( LSAME( UPLO, 'U' ) ) UP = .TRUE.
187: *
188: * Compute norm of op(A)*op2(C).
189: *
190: ANORM = 0.0D+0
191: IF ( UP ) THEN
192: DO I = 1, N
193: TMP = 0.0D+0
194: DO J = 1, I
195: TMP = TMP + CABS1( A( J, I ) * X( J ) )
196: END DO
197: DO J = I+1, N
198: TMP = TMP + CABS1( A( I, J ) * X( J ) )
199: END DO
200: RWORK( I ) = TMP
201: ANORM = MAX( ANORM, TMP )
202: END DO
203: ELSE
204: DO I = 1, N
205: TMP = 0.0D+0
206: DO J = 1, I
207: TMP = TMP + CABS1( A( I, J ) * X( J ) )
208: END DO
209: DO J = I+1, N
210: TMP = TMP + CABS1( A( J, I ) * X( J ) )
211: END DO
212: RWORK( I ) = TMP
213: ANORM = MAX( ANORM, TMP )
214: END DO
215: END IF
216: *
217: * Quick return if possible.
218: *
219: IF( N.EQ.0 ) THEN
220: ZLA_PORCOND_X = 1.0D+0
221: RETURN
222: ELSE IF( ANORM .EQ. 0.0D+0 ) THEN
223: RETURN
224: END IF
225: *
226: * Estimate the norm of inv(op(A)).
227: *
228: AINVNM = 0.0D+0
229: *
230: KASE = 0
231: 10 CONTINUE
232: CALL ZLACN2( N, WORK( N+1 ), WORK, AINVNM, KASE, ISAVE )
233: IF( KASE.NE.0 ) THEN
234: IF( KASE.EQ.2 ) THEN
235: *
236: * Multiply by R.
237: *
238: DO I = 1, N
239: WORK( I ) = WORK( I ) * RWORK( I )
240: END DO
241: *
242: IF ( UP ) THEN
243: CALL ZPOTRS( 'U', N, 1, AF, LDAF,
244: $ WORK, N, INFO )
245: ELSE
246: CALL ZPOTRS( 'L', N, 1, AF, LDAF,
247: $ WORK, N, INFO )
248: ENDIF
249: *
250: * Multiply by inv(X).
251: *
252: DO I = 1, N
253: WORK( I ) = WORK( I ) / X( I )
254: END DO
255: ELSE
256: *
257: * Multiply by inv(X**H).
258: *
259: DO I = 1, N
260: WORK( I ) = WORK( I ) / X( I )
261: END DO
262: *
263: IF ( UP ) THEN
264: CALL ZPOTRS( 'U', N, 1, AF, LDAF,
265: $ WORK, N, INFO )
266: ELSE
267: CALL ZPOTRS( 'L', N, 1, AF, LDAF,
268: $ WORK, N, INFO )
269: END IF
270: *
271: * Multiply by R.
272: *
273: DO I = 1, N
274: WORK( I ) = WORK( I ) * RWORK( I )
275: END DO
276: END IF
277: GO TO 10
278: END IF
279: *
280: * Compute the estimate of the reciprocal condition number.
281: *
282: IF( AINVNM .NE. 0.0D+0 )
283: $ ZLA_PORCOND_X = 1.0D+0 / AINVNM
284: *
285: RETURN
286: *
287: * End of ZLA_PORCOND_X
288: *
289: END
CVSweb interface <joel.bertrand@systella.fr>