LAPACK 3.12.0 LAPACK: Linear Algebra PACKage
Searching...
No Matches
dpocon.f
Go to the documentation of this file.
1*> \brief \b DPOCON
2*
3* =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6* http://www.netlib.org/lapack/explore-html/
7*
8*> \htmlonly
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dpocon.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dpocon.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dpocon.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18* Definition:
19* ===========
20*
21* SUBROUTINE DPOCON( UPLO, N, A, LDA, ANORM, RCOND, WORK, IWORK,
22* INFO )
23*
24* .. Scalar Arguments ..
25* CHARACTER UPLO
26* INTEGER INFO, LDA, N
27* DOUBLE PRECISION ANORM, RCOND
28* ..
29* .. Array Arguments ..
30* INTEGER IWORK( * )
31* DOUBLE PRECISION A( LDA, * ), WORK( * )
32* ..
33*
34*
35*> \par Purpose:
36* =============
37*>
38*> \verbatim
39*>
40*> DPOCON estimates the reciprocal of the condition number (in the
41*> 1-norm) of a real symmetric positive definite matrix using the
42*> Cholesky factorization A = U**T*U or A = L*L**T computed by DPOTRF.
43*>
44*> An estimate is obtained for norm(inv(A)), and the reciprocal of the
45*> condition number is computed as RCOND = 1 / (ANORM * norm(inv(A))).
46*> \endverbatim
47*
48* Arguments:
49* ==========
50*
51*> \param[in] UPLO
52*> \verbatim
53*> UPLO is CHARACTER*1
54*> = 'U': Upper triangle of A is stored;
55*> = 'L': Lower triangle of A is stored.
56*> \endverbatim
57*>
58*> \param[in] N
59*> \verbatim
60*> N is INTEGER
61*> The order of the matrix A. N >= 0.
62*> \endverbatim
63*>
64*> \param[in] A
65*> \verbatim
66*> A is DOUBLE PRECISION array, dimension (LDA,N)
67*> The triangular factor U or L from the Cholesky factorization
68*> A = U**T*U or A = L*L**T, as computed by DPOTRF.
69*> \endverbatim
70*>
71*> \param[in] LDA
72*> \verbatim
73*> LDA is INTEGER
74*> The leading dimension of the array A. LDA >= max(1,N).
75*> \endverbatim
76*>
77*> \param[in] ANORM
78*> \verbatim
79*> ANORM is DOUBLE PRECISION
80*> The 1-norm (or infinity-norm) of the symmetric matrix A.
81*> \endverbatim
82*>
83*> \param[out] RCOND
84*> \verbatim
85*> RCOND is DOUBLE PRECISION
86*> The reciprocal of the condition number of the matrix A,
87*> computed as RCOND = 1/(ANORM * AINVNM), where AINVNM is an
88*> estimate of the 1-norm of inv(A) computed in this routine.
89*> \endverbatim
90*>
91*> \param[out] WORK
92*> \verbatim
93*> WORK is DOUBLE PRECISION array, dimension (3*N)
94*> \endverbatim
95*>
96*> \param[out] IWORK
97*> \verbatim
98*> IWORK is INTEGER array, dimension (N)
99*> \endverbatim
100*>
101*> \param[out] INFO
102*> \verbatim
103*> INFO is INTEGER
104*> = 0: successful exit
105*> < 0: if INFO = -i, the i-th argument had an illegal value
106*> \endverbatim
107*
108* Authors:
109* ========
110*
111*> \author Univ. of Tennessee
112*> \author Univ. of California Berkeley
113*> \author Univ. of Colorado Denver
114*> \author NAG Ltd.
115*
116*> \ingroup pocon
117*
118* =====================================================================
119 SUBROUTINE dpocon( UPLO, N, A, LDA, ANORM, RCOND, WORK, IWORK,
120 \$ INFO )
121*
122* -- LAPACK computational routine --
123* -- LAPACK is a software package provided by Univ. of Tennessee, --
124* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
125*
126* .. Scalar Arguments ..
127 CHARACTER UPLO
128 INTEGER INFO, LDA, N
129 DOUBLE PRECISION ANORM, RCOND
130* ..
131* .. Array Arguments ..
132 INTEGER IWORK( * )
133 DOUBLE PRECISION A( LDA, * ), WORK( * )
134* ..
135*
136* =====================================================================
137*
138* .. Parameters ..
139 DOUBLE PRECISION ONE, ZERO
140 parameter( one = 1.0d+0, zero = 0.0d+0 )
141* ..
142* .. Local Scalars ..
143 LOGICAL UPPER
144 CHARACTER NORMIN
145 INTEGER IX, KASE
146 DOUBLE PRECISION AINVNM, SCALE, SCALEL, SCALEU, SMLNUM
147* ..
148* .. Local Arrays ..
149 INTEGER ISAVE( 3 )
150* ..
151* .. External Functions ..
152 LOGICAL LSAME
153 INTEGER IDAMAX
154 DOUBLE PRECISION DLAMCH
155 EXTERNAL lsame, idamax, dlamch
156* ..
157* .. External Subroutines ..
158 EXTERNAL dlacn2, dlatrs, drscl, xerbla
159* ..
160* .. Intrinsic Functions ..
161 INTRINSIC abs, max
162* ..
163* .. Executable Statements ..
164*
165* Test the input parameters.
166*
167 info = 0
168 upper = lsame( uplo, 'U' )
169 IF( .NOT.upper .AND. .NOT.lsame( uplo, 'L' ) ) THEN
170 info = -1
171 ELSE IF( n.LT.0 ) THEN
172 info = -2
173 ELSE IF( lda.LT.max( 1, n ) ) THEN
174 info = -4
175 ELSE IF( anorm.LT.zero ) THEN
176 info = -5
177 END IF
178 IF( info.NE.0 ) THEN
179 CALL xerbla( 'DPOCON', -info )
180 RETURN
181 END IF
182*
183* Quick return if possible
184*
185 rcond = zero
186 IF( n.EQ.0 ) THEN
187 rcond = one
188 RETURN
189 ELSE IF( anorm.EQ.zero ) THEN
190 RETURN
191 END IF
192*
193 smlnum = dlamch( 'Safe minimum' )
194*
195* Estimate the 1-norm of inv(A).
196*
197 kase = 0
198 normin = 'N'
199 10 CONTINUE
200 CALL dlacn2( n, work( n+1 ), work, iwork, ainvnm, kase, isave )
201 IF( kase.NE.0 ) THEN
202 IF( upper ) THEN
203*
204* Multiply by inv(U**T).
205*
206 CALL dlatrs( 'Upper', 'Transpose', 'Non-unit', normin, n, a,
207 \$ lda, work, scalel, work( 2*n+1 ), info )
208 normin = 'Y'
209*
210* Multiply by inv(U).
211*
212 CALL dlatrs( 'Upper', 'No transpose', 'Non-unit', normin, n,
213 \$ a, lda, work, scaleu, work( 2*n+1 ), info )
214 ELSE
215*
216* Multiply by inv(L).
217*
218 CALL dlatrs( 'Lower', 'No transpose', 'Non-unit', normin, n,
219 \$ a, lda, work, scalel, work( 2*n+1 ), info )
220 normin = 'Y'
221*
222* Multiply by inv(L**T).
223*
224 CALL dlatrs( 'Lower', 'Transpose', 'Non-unit', normin, n, a,
225 \$ lda, work, scaleu, work( 2*n+1 ), info )
226 END IF
227*
228* Multiply by 1/SCALE if doing so will not cause overflow.
229*
230 scale = scalel*scaleu
231 IF( scale.NE.one ) THEN
232 ix = idamax( n, work, 1 )
233 IF( scale.LT.abs( work( ix ) )*smlnum .OR. scale.EQ.zero )
234 \$ GO TO 20
235 CALL drscl( n, scale, work, 1 )
236 END IF
237 GO TO 10
238 END IF
239*
240* Compute the estimate of the reciprocal condition number.
241*
242 IF( ainvnm.NE.zero )
243 \$ rcond = ( one / ainvnm ) / anorm
244*
245 20 CONTINUE
246 RETURN
247*
248* End of DPOCON
249*
250 END
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine dlacn2(n, v, x, isgn, est, kase, isave)
DLACN2 estimates the 1-norm of a square matrix, using reverse communication for evaluating matrix-vec...
Definition dlacn2.f:136
subroutine dlatrs(uplo, trans, diag, normin, n, a, lda, x, scale, cnorm, info)
DLATRS solves a triangular system of equations with the scale factor set to prevent overflow.
Definition dlatrs.f:238
subroutine dpocon(uplo, n, a, lda, anorm, rcond, work, iwork, info)
DPOCON
Definition dpocon.f:121
subroutine drscl(n, sa, sx, incx)
DRSCL multiplies a vector by the reciprocal of a real scalar.
Definition drscl.f:84