LAPACK 3.12.1
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
dla_gercond.f
Go to the documentation of this file.
1*> \brief \b DLA_GERCOND estimates the Skeel condition number for a general matrix.
2*
3* =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6* http://www.netlib.org/lapack/explore-html/
7*
8*> Download DLA_GERCOND + dependencies
9*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dla_gercond.f">
10*> [TGZ]</a>
11*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dla_gercond.f">
12*> [ZIP]</a>
13*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dla_gercond.f">
14*> [TXT]</a>
15*
16* Definition:
17* ===========
18*
19* DOUBLE PRECISION FUNCTION DLA_GERCOND( TRANS, N, A, LDA, AF,
20* LDAF, IPIV, CMODE, C,
21* INFO, WORK, IWORK )
22*
23* .. Scalar Arguments ..
24* CHARACTER TRANS
25* INTEGER N, LDA, LDAF, INFO, CMODE
26* ..
27* .. Array Arguments ..
28* INTEGER IPIV( * ), IWORK( * )
29* DOUBLE PRECISION A( LDA, * ), AF( LDAF, * ), WORK( * ),
30* $ C( * )
31* ..
32*
33*
34*> \par Purpose:
35* =============
36*>
37*> \verbatim
38*>
39*> DLA_GERCOND estimates the Skeel condition number of op(A) * op2(C)
40*> where op2 is determined by CMODE as follows
41*> CMODE = 1 op2(C) = C
42*> CMODE = 0 op2(C) = I
43*> CMODE = -1 op2(C) = inv(C)
44*> The Skeel condition number cond(A) = norminf( |inv(A)||A| )
45*> is computed by computing scaling factors R such that
46*> diag(R)*A*op2(C) is row equilibrated and computing the standard
47*> infinity-norm condition number.
48*> \endverbatim
49*
50* Arguments:
51* ==========
52*
53*> \param[in] TRANS
54*> \verbatim
55*> TRANS is CHARACTER*1
56*> Specifies the form of the system of equations:
57*> = 'N': A * X = B (No transpose)
58*> = 'T': A**T * X = B (Transpose)
59*> = 'C': A**H * X = B (Conjugate Transpose = Transpose)
60*> \endverbatim
61*>
62*> \param[in] N
63*> \verbatim
64*> N is INTEGER
65*> The number of linear equations, i.e., the order of the
66*> matrix A. N >= 0.
67*> \endverbatim
68*>
69*> \param[in] A
70*> \verbatim
71*> A is DOUBLE PRECISION array, dimension (LDA,N)
72*> On entry, the N-by-N matrix A.
73*> \endverbatim
74*>
75*> \param[in] LDA
76*> \verbatim
77*> LDA is INTEGER
78*> The leading dimension of the array A. LDA >= max(1,N).
79*> \endverbatim
80*>
81*> \param[in] AF
82*> \verbatim
83*> AF is DOUBLE PRECISION array, dimension (LDAF,N)
84*> The factors L and U from the factorization
85*> A = P*L*U as computed by DGETRF.
86*> \endverbatim
87*>
88*> \param[in] LDAF
89*> \verbatim
90*> LDAF is INTEGER
91*> The leading dimension of the array AF. LDAF >= max(1,N).
92*> \endverbatim
93*>
94*> \param[in] IPIV
95*> \verbatim
96*> IPIV is INTEGER array, dimension (N)
97*> The pivot indices from the factorization A = P*L*U
98*> as computed by DGETRF; row i of the matrix was interchanged
99*> with row IPIV(i).
100*> \endverbatim
101*>
102*> \param[in] CMODE
103*> \verbatim
104*> CMODE is INTEGER
105*> Determines op2(C) in the formula op(A) * op2(C) as follows:
106*> CMODE = 1 op2(C) = C
107*> CMODE = 0 op2(C) = I
108*> CMODE = -1 op2(C) = inv(C)
109*> \endverbatim
110*>
111*> \param[in] C
112*> \verbatim
113*> C is DOUBLE PRECISION array, dimension (N)
114*> The vector C in the formula op(A) * op2(C).
115*> \endverbatim
116*>
117*> \param[out] INFO
118*> \verbatim
119*> INFO is INTEGER
120*> = 0: Successful exit.
121*> i > 0: The ith argument is invalid.
122*> \endverbatim
123*>
124*> \param[out] WORK
125*> \verbatim
126*> WORK is DOUBLE PRECISION array, dimension (3*N).
127*> Workspace.
128*> \endverbatim
129*>
130*> \param[out] IWORK
131*> \verbatim
132*> IWORK is INTEGER array, dimension (N).
133*> Workspace.
134*> \endverbatim
135*
136* Authors:
137* ========
138*
139*> \author Univ. of Tennessee
140*> \author Univ. of California Berkeley
141*> \author Univ. of Colorado Denver
142*> \author NAG Ltd.
143*
144*> \ingroup la_gercond
145*
146* =====================================================================
147 DOUBLE PRECISION FUNCTION dla_gercond( TRANS, N, A, LDA, AF,
148 $ LDAF, IPIV, CMODE, C,
149 $ INFO, WORK, IWORK )
150*
151* -- LAPACK computational routine --
152* -- LAPACK is a software package provided by Univ. of Tennessee, --
153* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
154*
155* .. Scalar Arguments ..
156 CHARACTER trans
157 INTEGER n, lda, ldaf, info, cmode
158* ..
159* .. Array Arguments ..
160 INTEGER ipiv( * ), iwork( * )
161 DOUBLE PRECISION a( lda, * ), af( ldaf, * ), work( * ),
162 $ c( * )
163* ..
164*
165* =====================================================================
166*
167* .. Local Scalars ..
168 LOGICAL notrans
169 INTEGER kase, i, j
170 DOUBLE PRECISION ainvnm, tmp
171* ..
172* .. Local Arrays ..
173 INTEGER isave( 3 )
174* ..
175* .. External Functions ..
176 LOGICAL lsame
177 EXTERNAL lsame
178* ..
179* .. External Subroutines ..
180 EXTERNAL dlacn2, dgetrs, xerbla
181* ..
182* .. Intrinsic Functions ..
183 INTRINSIC abs, max
184* ..
185* .. Executable Statements ..
186*
187 dla_gercond = 0.0d+0
188*
189 info = 0
190 notrans = lsame( trans, 'N' )
191 IF ( .NOT. notrans .AND. .NOT. lsame(trans, 'T')
192 $ .AND. .NOT. lsame(trans, 'C') ) THEN
193 info = -1
194 ELSE IF( n.LT.0 ) THEN
195 info = -2
196 ELSE IF( lda.LT.max( 1, n ) ) THEN
197 info = -4
198 ELSE IF( ldaf.LT.max( 1, n ) ) THEN
199 info = -6
200 END IF
201 IF( info.NE.0 ) THEN
202 CALL xerbla( 'DLA_GERCOND', -info )
203 RETURN
204 END IF
205 IF( n.EQ.0 ) THEN
206 dla_gercond = 1.0d+0
207 RETURN
208 END IF
209*
210* Compute the equilibration matrix R such that
211* inv(R)*A*C has unit 1-norm.
212*
213 IF (notrans) THEN
214 DO i = 1, n
215 tmp = 0.0d+0
216 IF ( cmode .EQ. 1 ) THEN
217 DO j = 1, n
218 tmp = tmp + abs( a( i, j ) * c( j ) )
219 END DO
220 ELSE IF ( cmode .EQ. 0 ) THEN
221 DO j = 1, n
222 tmp = tmp + abs( a( i, j ) )
223 END DO
224 ELSE
225 DO j = 1, n
226 tmp = tmp + abs( a( i, j ) / c( j ) )
227 END DO
228 END IF
229 work( 2*n+i ) = tmp
230 END DO
231 ELSE
232 DO i = 1, n
233 tmp = 0.0d+0
234 IF ( cmode .EQ. 1 ) THEN
235 DO j = 1, n
236 tmp = tmp + abs( a( j, i ) * c( j ) )
237 END DO
238 ELSE IF ( cmode .EQ. 0 ) THEN
239 DO j = 1, n
240 tmp = tmp + abs( a( j, i ) )
241 END DO
242 ELSE
243 DO j = 1, n
244 tmp = tmp + abs( a( j, i ) / c( j ) )
245 END DO
246 END IF
247 work( 2*n+i ) = tmp
248 END DO
249 END IF
250*
251* Estimate the norm of inv(op(A)).
252*
253 ainvnm = 0.0d+0
254
255 kase = 0
256 10 CONTINUE
257 CALL dlacn2( n, work( n+1 ), work, iwork, ainvnm, kase, isave )
258 IF( kase.NE.0 ) THEN
259 IF( kase.EQ.2 ) THEN
260*
261* Multiply by R.
262*
263 DO i = 1, n
264 work(i) = work(i) * work(2*n+i)
265 END DO
266
267 IF (notrans) THEN
268 CALL dgetrs( 'No transpose', n, 1, af, ldaf, ipiv,
269 $ work, n, info )
270 ELSE
271 CALL dgetrs( 'Transpose', n, 1, af, ldaf, ipiv,
272 $ work, n, info )
273 END IF
274*
275* Multiply by inv(C).
276*
277 IF ( cmode .EQ. 1 ) THEN
278 DO i = 1, n
279 work( i ) = work( i ) / c( i )
280 END DO
281 ELSE IF ( cmode .EQ. -1 ) THEN
282 DO i = 1, n
283 work( i ) = work( i ) * c( i )
284 END DO
285 END IF
286 ELSE
287*
288* Multiply by inv(C**T).
289*
290 IF ( cmode .EQ. 1 ) THEN
291 DO i = 1, n
292 work( i ) = work( i ) / c( i )
293 END DO
294 ELSE IF ( cmode .EQ. -1 ) THEN
295 DO i = 1, n
296 work( i ) = work( i ) * c( i )
297 END DO
298 END IF
299
300 IF (notrans) THEN
301 CALL dgetrs( 'Transpose', n, 1, af, ldaf, ipiv,
302 $ work, n, info )
303 ELSE
304 CALL dgetrs( 'No transpose', n, 1, af, ldaf, ipiv,
305 $ work, n, info )
306 END IF
307*
308* Multiply by R.
309*
310 DO i = 1, n
311 work( i ) = work( i ) * work( 2*n+i )
312 END DO
313 END IF
314 GO TO 10
315 END IF
316*
317* Compute the estimate of the reciprocal condition number.
318*
319 IF( ainvnm .NE. 0.0d+0 )
320 $ dla_gercond = ( 1.0d+0 / ainvnm )
321*
322 RETURN
323*
324* End of DLA_GERCOND
325*
326 END
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine dgetrs(trans, n, nrhs, a, lda, ipiv, b, ldb, info)
DGETRS
Definition dgetrs.f:119
double precision function dla_gercond(trans, n, a, lda, af, ldaf, ipiv, cmode, c, info, work, iwork)
DLA_GERCOND estimates the Skeel condition number for a general matrix.
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:134
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48