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