LAPACK 3.12.1
LAPACK: Linear Algebra PACKage
•All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Modules Pages

◆ sla_porcond()

real function sla_porcond ( character uplo,
integer n,
real, dimension( lda, * ) a,
integer lda,
real, dimension( ldaf, * ) af,
integer ldaf,
integer cmode,
real, dimension( * ) c,
integer info,
real, dimension( * ) work,
integer, dimension( * ) iwork )

SLA_PORCOND estimates the Skeel condition number for a symmetric positive-definite matrix.

Download SLA_PORCOND + dependencies [TGZ] [ZIP] [TXT]

Purpose:
!>
!>    SLA_PORCOND Estimates the Skeel condition number of  op(A) * op2(C)
!>    where op2 is determined by CMODE as follows
!>    CMODE =  1    op2(C) = C
!>    CMODE =  0    op2(C) = I
!>    CMODE = -1    op2(C) = inv(C)
!>    The Skeel condition number  cond(A) = norminf( |inv(A)||A| )
!>    is computed by computing scaling factors R such that
!>    diag(R)*A*op2(C) is row equilibrated and computing the standard
!>    infinity-norm condition number.
!>
Parameters
[in]UPLO
!>          UPLO is CHARACTER*1
!>       = 'U':  Upper triangle of A is stored;
!>       = 'L':  Lower triangle of A is stored.
!>
[in]N
!>          N is INTEGER
!>     The number of linear equations, i.e., the order of the
!>     matrix A.  N >= 0.
!>
[in]A
!>          A is REAL array, dimension (LDA,N)
!>     On entry, the N-by-N matrix A.
!>
[in]LDA
!>          LDA is INTEGER
!>     The leading dimension of the array A.  LDA >= max(1,N).
!>
[in]AF
!>          AF is REAL array, dimension (LDAF,N)
!>     The triangular factor U or L from the Cholesky factorization
!>     A = U**T*U or A = L*L**T, as computed by SPOTRF.
!>
[in]LDAF
!>          LDAF is INTEGER
!>     The leading dimension of the array AF.  LDAF >= max(1,N).
!>
[in]CMODE
!>          CMODE is INTEGER
!>     Determines op2(C) in the formula op(A) * op2(C) as follows:
!>     CMODE =  1    op2(C) = C
!>     CMODE =  0    op2(C) = I
!>     CMODE = -1    op2(C) = inv(C)
!>
[in]C
!>          C is REAL array, dimension (N)
!>     The vector C in the formula op(A) * op2(C).
!>
[out]INFO
!>          INFO is INTEGER
!>       = 0:  Successful exit.
!>     i > 0:  The ith argument is invalid.
!>
[out]WORK
!>          WORK is REAL array, dimension (3*N).
!>     Workspace.
!>
[out]IWORK
!>          IWORK is INTEGER array, dimension (N).
!>     Workspace.
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 136 of file sla_porcond.f.

138*
139* -- LAPACK computational routine --
140* -- LAPACK is a software package provided by Univ. of Tennessee, --
141* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
142*
143* .. Scalar Arguments ..
144 CHARACTER UPLO
145 INTEGER N, LDA, LDAF, INFO, CMODE
146 REAL A( LDA, * ), AF( LDAF, * ), WORK( * ),
147 $ C( * )
148* ..
149* .. Array Arguments ..
150 INTEGER IWORK( * )
151* ..
152*
153* =====================================================================
154*
155* .. Local Scalars ..
156 INTEGER KASE, I, J
157 REAL AINVNM, TMP
158 LOGICAL UP
159* ..
160* .. Array Arguments ..
161 INTEGER ISAVE( 3 )
162* ..
163* .. External Functions ..
164 LOGICAL LSAME
165 EXTERNAL lsame
166* ..
167* .. External Subroutines ..
168 EXTERNAL slacn2, spotrs, xerbla
169* ..
170* .. Intrinsic Functions ..
171 INTRINSIC abs, max
172* ..
173* .. Executable Statements ..
174*
175 sla_porcond = 0.0
176*
177 info = 0
178 IF( n.LT.0 ) THEN
179 info = -2
180 END IF
181 IF( info.NE.0 ) THEN
182 CALL xerbla( 'SLA_PORCOND', -info )
183 RETURN
184 END IF
185
186 IF( n.EQ.0 ) THEN
187 sla_porcond = 1.0
188 RETURN
189 END IF
190 up = .false.
191 IF ( lsame( uplo, 'U' ) ) up = .true.
192*
193* Compute the equilibration matrix R such that
194* inv(R)*A*C has unit 1-norm.
195*
196 IF ( up ) THEN
197 DO i = 1, n
198 tmp = 0.0
199 IF ( cmode .EQ. 1 ) THEN
200 DO j = 1, i
201 tmp = tmp + abs( a( j, i ) * c( j ) )
202 END DO
203 DO j = i+1, n
204 tmp = tmp + abs( a( i, j ) * c( j ) )
205 END DO
206 ELSE IF ( cmode .EQ. 0 ) THEN
207 DO j = 1, i
208 tmp = tmp + abs( a( j, i ) )
209 END DO
210 DO j = i+1, n
211 tmp = tmp + abs( a( i, j ) )
212 END DO
213 ELSE
214 DO j = 1, i
215 tmp = tmp + abs( a( j ,i ) / c( j ) )
216 END DO
217 DO j = i+1, n
218 tmp = tmp + abs( a( i, j ) / c( j ) )
219 END DO
220 END IF
221 work( 2*n+i ) = tmp
222 END DO
223 ELSE
224 DO i = 1, n
225 tmp = 0.0
226 IF ( cmode .EQ. 1 ) THEN
227 DO j = 1, i
228 tmp = tmp + abs( a( i, j ) * c( j ) )
229 END DO
230 DO j = i+1, n
231 tmp = tmp + abs( a( j, i ) * c( j ) )
232 END DO
233 ELSE IF ( cmode .EQ. 0 ) THEN
234 DO j = 1, i
235 tmp = tmp + abs( a( i, j ) )
236 END DO
237 DO j = i+1, n
238 tmp = tmp + abs( a( j, i ) )
239 END DO
240 ELSE
241 DO j = 1, i
242 tmp = tmp + abs( a( i, j ) / c( j ) )
243 END DO
244 DO j = i+1, n
245 tmp = tmp + abs( a( j, i ) / c( j ) )
246 END DO
247 END IF
248 work( 2*n+i ) = tmp
249 END DO
250 ENDIF
251*
252* Estimate the norm of inv(op(A)).
253*
254 ainvnm = 0.0
255
256 kase = 0
257 10 CONTINUE
258 CALL slacn2( n, work( n+1 ), work, iwork, ainvnm, kase, isave )
259 IF( kase.NE.0 ) THEN
260 IF( kase.EQ.2 ) THEN
261*
262* Multiply by R.
263*
264 DO i = 1, n
265 work( i ) = work( i ) * work( 2*n+i )
266 END DO
267
268 IF (up) THEN
269 CALL spotrs( 'Upper', n, 1, af, ldaf, work, n, info )
270 ELSE
271 CALL spotrs( 'Lower', n, 1, af, ldaf, work, n, info )
272 ENDIF
273*
274* Multiply by inv(C).
275*
276 IF ( cmode .EQ. 1 ) THEN
277 DO i = 1, n
278 work( i ) = work( i ) / c( i )
279 END DO
280 ELSE IF ( cmode .EQ. -1 ) THEN
281 DO i = 1, n
282 work( i ) = work( i ) * c( i )
283 END DO
284 END IF
285 ELSE
286*
287* Multiply by inv(C**T).
288*
289 IF ( cmode .EQ. 1 ) THEN
290 DO i = 1, n
291 work( i ) = work( i ) / c( i )
292 END DO
293 ELSE IF ( cmode .EQ. -1 ) THEN
294 DO i = 1, n
295 work( i ) = work( i ) * c( i )
296 END DO
297 END IF
298
299 IF ( up ) THEN
300 CALL spotrs( 'Upper', n, 1, af, ldaf, work, n, info )
301 ELSE
302 CALL spotrs( 'Lower', n, 1, af, ldaf, work, n, info )
303 ENDIF
304*
305* Multiply by R.
306*
307 DO i = 1, n
308 work( i ) = work( i ) * work( 2*n+i )
309 END DO
310 END IF
311 GO TO 10
312 END IF
313*
314* Compute the estimate of the reciprocal condition number.
315*
316 IF( ainvnm .NE. 0.0 )
317 $ sla_porcond = ( 1.0 / ainvnm )
318*
319 RETURN
320*
321* End of SLA_PORCOND
322*
xerbla
subroutine xerbla(srname, info)
Definition cblat2.f:3285
sla_porcond
real function sla_porcond(uplo, n, a, lda, af, ldaf, cmode, c, info, work, iwork)
SLA_PORCOND estimates the Skeel condition number for a symmetric positive-definite matrix.
Definition sla_porcond.f:138
slacn2
subroutine slacn2(n, v, x, isgn, est, kase, isave)
SLACN2 estimates the 1-norm of a square matrix, using reverse communication for evaluating matrix-vec...
Definition slacn2.f:134
lsame
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
spotrs
subroutine spotrs(uplo, n, nrhs, a, lda, b, ldb, info)
SPOTRS
Definition spotrs.f:108
Here is the call graph for this function:
Here is the caller graph for this function: