Annotation of rpl/lapack/lapack/zhpcon.f, revision 1.5
1.1 bertrand 1: SUBROUTINE ZHPCON( UPLO, N, AP, IPIV, ANORM, RCOND, WORK, INFO )
2: *
3: * -- LAPACK routine (version 3.2) --
4: * -- LAPACK is a software package provided by Univ. of Tennessee, --
5: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
6: * November 2006
7: *
8: * Modified to call ZLACN2 in place of ZLACON, 10 Feb 03, SJH.
9: *
10: * .. Scalar Arguments ..
11: CHARACTER UPLO
12: INTEGER INFO, N
13: DOUBLE PRECISION ANORM, RCOND
14: * ..
15: * .. Array Arguments ..
16: INTEGER IPIV( * )
17: COMPLEX*16 AP( * ), WORK( * )
18: * ..
19: *
20: * Purpose
21: * =======
22: *
23: * ZHPCON estimates the reciprocal of the condition number of a complex
24: * Hermitian packed matrix A using the factorization A = U*D*U**H or
25: * A = L*D*L**H computed by ZHPTRF.
26: *
27: * An estimate is obtained for norm(inv(A)), and the reciprocal of the
28: * condition number is computed as RCOND = 1 / (ANORM * norm(inv(A))).
29: *
30: * Arguments
31: * =========
32: *
33: * UPLO (input) CHARACTER*1
34: * Specifies whether the details of the factorization are stored
35: * as an upper or lower triangular matrix.
36: * = 'U': Upper triangular, form is A = U*D*U**H;
37: * = 'L': Lower triangular, form is A = L*D*L**H.
38: *
39: * N (input) INTEGER
40: * The order of the matrix A. N >= 0.
41: *
42: * AP (input) COMPLEX*16 array, dimension (N*(N+1)/2)
43: * The block diagonal matrix D and the multipliers used to
44: * obtain the factor U or L as computed by ZHPTRF, stored as a
45: * packed triangular matrix.
46: *
47: * IPIV (input) INTEGER array, dimension (N)
48: * Details of the interchanges and the block structure of D
49: * as determined by ZHPTRF.
50: *
51: * ANORM (input) DOUBLE PRECISION
52: * The 1-norm of the original matrix A.
53: *
54: * RCOND (output) DOUBLE PRECISION
55: * The reciprocal of the condition number of the matrix A,
56: * computed as RCOND = 1/(ANORM * AINVNM), where AINVNM is an
57: * estimate of the 1-norm of inv(A) computed in this routine.
58: *
59: * WORK (workspace) COMPLEX*16 array, dimension (2*N)
60: *
61: * INFO (output) INTEGER
62: * = 0: successful exit
63: * < 0: if INFO = -i, the i-th argument had an illegal value
64: *
65: * =====================================================================
66: *
67: * .. Parameters ..
68: DOUBLE PRECISION ONE, ZERO
69: PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 )
70: * ..
71: * .. Local Scalars ..
72: LOGICAL UPPER
73: INTEGER I, IP, KASE
74: DOUBLE PRECISION AINVNM
75: * ..
76: * .. Local Arrays ..
77: INTEGER ISAVE( 3 )
78: * ..
79: * .. External Functions ..
80: LOGICAL LSAME
81: EXTERNAL LSAME
82: * ..
83: * .. External Subroutines ..
84: EXTERNAL XERBLA, ZHPTRS, ZLACN2
85: * ..
86: * .. Executable Statements ..
87: *
88: * Test the input parameters.
89: *
90: INFO = 0
91: UPPER = LSAME( UPLO, 'U' )
92: IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
93: INFO = -1
94: ELSE IF( N.LT.0 ) THEN
95: INFO = -2
96: ELSE IF( ANORM.LT.ZERO ) THEN
97: INFO = -5
98: END IF
99: IF( INFO.NE.0 ) THEN
100: CALL XERBLA( 'ZHPCON', -INFO )
101: RETURN
102: END IF
103: *
104: * Quick return if possible
105: *
106: RCOND = ZERO
107: IF( N.EQ.0 ) THEN
108: RCOND = ONE
109: RETURN
110: ELSE IF( ANORM.LE.ZERO ) THEN
111: RETURN
112: END IF
113: *
114: * Check that the diagonal matrix D is nonsingular.
115: *
116: IF( UPPER ) THEN
117: *
118: * Upper triangular storage: examine D from bottom to top
119: *
120: IP = N*( N+1 ) / 2
121: DO 10 I = N, 1, -1
122: IF( IPIV( I ).GT.0 .AND. AP( IP ).EQ.ZERO )
123: $ RETURN
124: IP = IP - I
125: 10 CONTINUE
126: ELSE
127: *
128: * Lower triangular storage: examine D from top to bottom.
129: *
130: IP = 1
131: DO 20 I = 1, N
132: IF( IPIV( I ).GT.0 .AND. AP( IP ).EQ.ZERO )
133: $ RETURN
134: IP = IP + N - I + 1
135: 20 CONTINUE
136: END IF
137: *
138: * Estimate the 1-norm of the inverse.
139: *
140: KASE = 0
141: 30 CONTINUE
142: CALL ZLACN2( N, WORK( N+1 ), WORK, AINVNM, KASE, ISAVE )
143: IF( KASE.NE.0 ) THEN
144: *
145: * Multiply by inv(L*D*L') or inv(U*D*U').
146: *
147: CALL ZHPTRS( UPLO, N, 1, AP, IPIV, WORK, N, INFO )
148: GO TO 30
149: END IF
150: *
151: * Compute the estimate of the reciprocal condition number.
152: *
153: IF( AINVNM.NE.ZERO )
154: $ RCOND = ( ONE / AINVNM ) / ANORM
155: *
156: RETURN
157: *
158: * End of ZHPCON
159: *
160: END
CVSweb interface <joel.bertrand@systella.fr>
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to zhpcon.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / lapack